Global CO2 fluxes estimated from GOSAT retrievals of total column CO2

We present one of the first estimates of the global distribution of CO2 surface fluxes using total column CO 2 measurements retrieved by the SRON-KIT RemoTeC algorithm from the Greenhouse gases Observing SATellite (GOSAT). We derive optimized fluxes from June 2009 to December 2010. We estimate fluxes from surface CO 2 measurements to use as baselines for comparing GOSAT dataderived fluxes. Assimilating only GOSAT data, we can reproduce the observed CO 2 time series at surface and TCCON sites in the tropics and the northern extra-tropics. In contrast, in the southern extra-tropics GOSAT X CO2 leads to enhanced seasonal cycle amplitudes compared to independent measurements, and we identify it as the result of a land– sea bias in our GOSAT X CO2 retrievals. A bias correction in the form of a global offset between GOSAT land and sea pixels in a joint inversion of satellite and surface measurements of CO2 yields plausible global flux estimates which are more tightly constrained than in an inversion using surface CO2 data alone. We show that assimilating the biascorrected GOSAT data on top of surface CO 2 data (a) reduces the estimated global land sink of CO 2, and (b) shifts the terrestrial net uptake of carbon from the tropics to the extratropics. It is concluded that while GOSAT total column CO 2 provide useful constraints for source–sink inversions, small spatiotemporal biases – beyond what can be detected using current validation techniques – have serious consequences for optimized fluxes, even aggregated over continental scales.


Introduction
The traditional top-down approach for quantifying sources and sinks of CO 2 is an atmospheric inversion of CO 2 concentrations measured at the earth's surface (Rödenbeck et al., 2003;Peters et al., 2007;Mueller et al., 2008;Chevallier et al., 2010;Gurney et al., 2008;Gurney, 2004). The measurements are performed by a network of trace gas monitoring stations established by various national agencies (Tans et al., 1990;Beardsmode and Pearman, 1987;Maki et al., 2011;Biraud et al., 2013). For any given inversion framework, the uncertainty bounds on regional CO 2 source-sink estimates depend heavily on the density of surface stations. The sparseness and spatial inhomogeneity of the existing surface network have limited our ability to understand the quantity and spatiotemporal distribution of CO 2 sources and sinks (Scholes et al., 2009). The difficulty of establishing surface stations in many areas that are interesting from a carbon cycle point of view -such as the Amazonian rain forest and the Arctic tundra -has prompted the development of spacebased CO 2 sensors. Currently operational instruments such as the Infrared Atmospheric Sounding Interferometer (IASI, Crevoisier et al., 2009), the Operational Vertical Sounder (TOVS, Chédin et al., 2003), the Tropospheric Emissions Spectrometer (TES, Kulawik et al., 2010), and the Atmospheric Infrared Sounder (AIRS, Chahine et al., 2006), although sensitive to atmospheric CO 2 , were not designed to study surface sources and sinks of CO 2 , and lack sensitivity to near-surface CO 2 . The Scanning Imaging Absorption Spectrometer for Atmospheric Cartography (SCIAMACHY) on board the ENVISAT satellite was the first space-based instrument to be sensitive to CO 2 in the lower troposphere (Reuter et al., 2011). More recent missions such as the Greenhouse gases Observing SATellite (GOSAT, Hamazaki et al., 2004) and planned missions such as the Orbiting Carbon Observatory-2 (OCO2, Boesch et al., 2011) and the Active Sensing of CO 2 Emissions over Nights, Days, and Seasons (ASCENDS, Dobbs et al., 2007) have been and are being developed specifically to resolve surface sources and sinks of CO 2 .
Several observation system simulation experiments have explored the added benefit of assimilating satellite measurements in atmospheric inversions of CO 2 (Rayner and O'Brien, 2001;Park and Prather, 2001;Houweling et al., 2004;Chevallier et al., 2007;Miller et al., 2007;Hungershoefer et al., 2010). These studies have shown that although space-based measurements of CO 2 are not as accurate as surface layer measurements, the increased spatial coverage can provide information not available from the sparse surface network (Buchwitz et al., 2007). Whether this information can be harvested for source-sink inversions, however, depends on the specifics of the satellite instrument. Variations in CO 2 concentrations caused by surface fluxes are damped by vertical transport, so satellite measurements that are sensitive to CO 2 at lower altitudes are better at constraining fluxes than measurements of CO 2 at higher altitudes. For example, Chevallier et al. showed that assimilating TOVS CO 2 observations in an inversion yielded unrealistic surface fluxes (Chevallier et al., 2005), and their attempt to assimilate AIRS radiances was outperformed by a standard surface data inversion (Chevallier et al., 2009). The poor performance of AIRS and TOVS in terms of constraining CO 2 surface fluxes can be traced to the fact that both these instruments are most sensitive to CO 2 in the upper troposphere (∼ 150 hPa). More recently, Nassar et al. (2011) tried assimilating CO 2 measurements from the TES instrument, which, with a spectral coverage of 3.2-15.4 µm, is sensitive to the mid-troposphere (∼ 550 hPa). They showed that assimilating TES CO 2 in addition to surface CO 2 measurements improves constraints on surface fluxes, especially over regions absent from the surface network, such as the tropical forests of South America (Nassar et al., 2011). It follows that a satellite instrument sensitive to the lower troposphere would be even more useful for CO 2 source-sink inversions.
The Thermal And Near infrared Sensor for carbon Observation (TANSO) instrument on board the GOSAT satellite (Kuze et al., 2009) acquires CO 2 spectra in the 1.6 µm and 2.0 µm bands. It is therefore sensitive to the CO 2 concentration in the mid-to lower troposphere, and is thus sensitive to surface fluxes of CO 2 . GOSAT has a polar orbit with a local overpass time of ∼ 13:00 and a three-day repeat cycle. TANSO measures the intensity of reflected sunlight at the two CO 2 absorption bands mentioned before, and is therefore sensitive to the total column CO 2 between the surface and the top of the atmosphere within its footprint area of ∼ 80 km 2 . We use the SRON-KIT RemoTeC retrieval algorithm to translate level 1 radiances to level 2 columnaveraged dry air CO 2 mole fractions, hereafter referred to as X CO 2 (Butz et al., 2011). RemoTeC is a "full physics" algorithm that simultaneously retrieves X CO 2 , X CH 4 and aerosol parameters needed to correct the optical path for the impact of scattering. The algorithm was evaluated by comparing the retrieved X CO 2 and X CH 4 with surface-based measurements of X CO 2 and X CH 4 from 12 stations of the Total Carbon Column Observing Network (TCCON, Wunch et al., 2011). Compared to TCCON X CO 2 , RemoTeC X CO 2 has a singleshot precision of 2.5 ppm and a bias of −0.36 % averaged over all stations .
In this manuscript, we investigate the value of assimilating GOSAT X CO 2 in a source-sink inversion of CO 2 using RemoTeC retrievals of X CO 2 . We investigate whether there are obvious biases between RemoTeC X CO 2 and surface CO 2 measurements, and devise strategies to correct for those biases to jointly assimilate X CO 2 and surface CO 2 data. We check if the addition of X CO 2 constraints on top of surface CO 2 data reduces the uncertainty on estimated fluxes from an atmospheric inversion, especially over areas with scant surface data. Most importantly, we look for portions of the carbon cycle about which our knowledge changes quantitatively when we jointly assimilate X CO 2 and surface CO 2 , compared to a surface-only inversion.
This idea of jointly assimilating observations from satellites and surface networks has been found to be more beneficial than assimilating either type of observation individually (Hungershoefer et al., 2010;Chevallier et al., 2009). Surface measurements of CO 2 are in general more accurate than satellite measurements, so a given number of surface observations provide a tighter constraint on surface fluxes than the same number of satellite observations. Satellite measurements, on the other hand, have a wider coverage, and therefore provide constraints on CO 2 in areas devoid of surface stations. Even in areas with a high density of surface stations, surface CO 2 and X CO 2 measurements can impose different constraints on surface fluxes, resulting in, for example, different amplitudes for the estimated seasonal cycle in the net ecosystem exchange (Yang et al., 2007;Basu et al., 2011;Keppel-Aleks et al., 2012).
A joint assimilation of surface CO 2 and satellite X CO 2 measurements can potentially compensate for spatiotemporal biases in the satellite instrument using ground data from background stations. Since most space-based CO 2 sensors are one of a kind, each comes with its own bias characteristics, which can vary seasonally, and with latitude, and can drift during the lifetime of the instrument. These biases can be partly compensated by assimilating background surface stations -such as Mauna Loa and South Pole -which provide accurate information on large-scale spatial gradients and long-term temporal trends of CO 2 .
The GOSAT X CO 2 data record starts from April 2009. In this work, we assimilate both surface and GOSAT data between 1 June 2009 and 1 January 2011, to optimize the surface flux of CO 2 between 1 June 2009 and 1 December 2010. The structure of this manuscript is as follows. In Sect. 2, we describe the satellite and surface data assimilated in our inversions. In Sect. 2.1 we describe the selection procedure for surface measurements and the temporal averaging employed. In Sect. 2.2 we describe the X CO 2 retrievals used, and in particular how errors are ascribed to individual soundings, which is crucial for striking the right balance between satellite and surface data. In Sect. 3 we describe the TM5 4DVAR atmospheric inversion system, with particular focus on the prior fluxes used in Sect. 3.1.1, the initial atmospheric CO 2 field in Sect. 3.1.3, and the estimation of the posterior covariance matrix in Sect. 3.2. Section 3.3 describes how we validate our optimized fluxes by comparing the resultant atmospheric CO 2 signal with aircraft measurements from the CONTRAIL campaign and X CO 2 measurements from the TCCON network. Section 4 contains the results of our atmospheric inversions, with Sect. 4.2 devoted to validation. Sensitivity runs to test the robustness of our optimization approach are described in the Supplement Sect. S1. Section 4.3 is where we present the estimated fluxes themselves, including large-scale aggregates and seasonal time series. Finally, we discuss the implications of our findings in Sect. 5.

Assimilated data
We assimilate both boundary layer CO 2 measurements, i.e., from surface flasks and in situ measurements, and GOSAT X CO 2 . In order to derive optimized fluxes from 1 June 2009 to 1 December 2010, we added three months of spin-up and one month of spin-down time, resulting in a model run from 1 March 2009 to 1 January 2011. However, we did not put any observational constraints during the spin-up time, and thus our assimilated data -both surface and satellite -span the period from 1 June 2009 to 1 January 2011. The two data streams are described in detail below.

Point samples
We assimilate 16 887 CO 2 surface layer observationsreferred to as "point samples" since each one represents the CO 2 concentration at a point in space and time -chosen from the set of observations assimilated in CarbonTracker 2011 (http://carbontracker.noaa.gov, ob- The marker size represents the amount of observational data from each location assimilated in our inversions. For example, although there is a high density of location markers over the tropical Pacific, there are very few data points from each of those locations, since the data there come from one-time cruises. servations retrieved from ftp://aftp.cmdl.noaa.gov/products/ carbontracker/co2/CT2011/observations/). The observations come from 132 distinct locations, which are shown in Fig. 1. The samples include data from nine tall towers (all within the continental US), thirteen in situ monitoring stations (seven from Environment Canada, two from the National Center for Atmospheric Research, and four from the National Oceanic and Atmospheric Administration), and sixty-three surface flask sites operated by various agencies. For tall towers having multiple intake levels, only the concentrations sampled from the highest levels are assimilated , such as the 396 m (aboveground) level of the WLEF tower near Park Falls, Wisconsin. Hourly CO 2 records were available from several CO 2 monitoring sites. However, assimilating the full hourly CO 2 record leads to model biases, especially at night (Patra et al., 2008). Following the prescription of Peters et al. (2010), we create at most one observation per station per day, representative of the largest possible footprint for that station. For lowaltitude stations, observations with the largest footprints are those collected when the daytime planetary boundary layer (PBL) is well developed and at its highest, which is during the mid to late afternoon. For mountaintop sites on the other hand, the largest footprint corresponds to the time when the observations sample free tropospheric air and the nighttime PBL is at its lowest, i.e., after midnight. Hence, depending on the station location, we assimilate either late afternoon or late night CO 2 data. Data collected are averaged within a four-hour time window to get rid of high frequency fluctuations that our coarse transport model cannot possibly resolve. During the data assimilation phase, the modeled CO 2 concentration at a site is also averaged over the same four-hour period before being compared to observations from that site.

8698
S. Basu et al.: Global CO 2 flux estimation using GOSAT X CO 2 2.2 Satellite measurements X CO 2 data used in this work were retrieved from GOSAT soundings using the RemoTeC algorithm. RemoTeC is being jointly developed at the Netherlands Institute for Space Research (SRON) and the Karlsruhe Institute of Technology (KIT). Its performance has been extensively evaluated through retrieval simulations, in particular with respect to its capability to account for lightpath modification by particle scattering in the atmosphere (Butz et al., 2010(Butz et al., , 2009. Initial validation against ground-based observations by TCCON has shown very promising data quality with respect to precision and accuracy of data. The performance of RemoTeC in comparison with other algorithms for GOSAT X CO 2 retrieval can be found in Oshchepkov et al. (2013).
Out of the 262 355 X CO 2 soundings retrieved by RemoTeC within our observation window, we assimilated 77 769 that passed RemoTeC quality checks. The assimilated X CO 2 observations, shown in Fig. 2 to give an idea of the coverage, include both land and ocean soundings, where the latter were acquired in GOSAT's ocean glint observation mode. Re-moTeC quality checks rely on a series of criteria described in Butz et al. (2011). GOSAT observations are screened for instrument malfunction, suitable observation geometry, signalto-noise, (cirrus) cloud contamination, and surface elevation variability of the observed scene before they are submitted to the retrieval. Data recorded in the "medium gain" mode of TANSO-FTS -typically collected over bright surfaces such as deserts -are discarded since a scan speed instability of the FTS affects data quality. Together with X CO 2 (and X CH 4 ), RemoTeC delivers a range of retrieval diagnostics that allow for a posteriori quality assessment, among them standard criteria such as the posterior χ 2 and degrees of freedom for signal. The retrieved particle scattering properties of the atmosphere hint at the difficulty of the observed scene in terms of lightpath modification through aerosol and cirrus particles. The retrieval diagnostics, along with the retrieved particle scattering properties, were used to select the 77 769 scenes mentioned above which were assimilated. Before data assimilation, a bias correction based on coincident TCCON soundings was applied to the 77 769 X CO 2 measurements, as detailed by Guerlet et al. (2013).
The 77 769 assimilated soundings do not necessarily contribute 77 769 independent observations. Errors in different X CO 2 retrievals could be correlated, for example, due to errors in ancillary data -such as meteorological fields and spectral line shapes -used in retrievals (see, e.g., Fig. 2 of Chevallier et al., 2005). Since observed X CO 2 values are assimilated by comparing them with simulated X CO 2 from an atmospheric transport model, the transport model can also generate correlations between model-observation mismatches (Kaminski et al., 2001). While accounting properly for these correlations in a Bayesian inversion is crucial, in reality the correlations are hard to estimate, and difficult to implement in an inversion framework. In lieu of a full obser-  The averaging is only for visual clarity; the data assimilation considered each observation individually. The color of a grid box represents the average X CO 2 within that grid box, in mole fraction of dry air expressed as parts per million. vation error covariance matrix, which is formally n obs × n obs where n obs ∼ 10 5 is the number of observations assimilated, data assimilation systems typically use any of a variety of techniques such as binning, data thinning and error inflation to construct a diagonal approximation of the full matrix, which is easier to implement. Of these approximations, error inflation was shown to perform best in the context of assimilating satellite X CO 2 for estimating CO 2 sources and sinks (Chevallier, 2007). They showed that multiplying all observation errors by 2 yielded optimized surface fluxes closest to the exact solution using the full error covariance matrix. In our inversion system, we also inflate X CO 2 errors to account for correlated observation errors. However, unlike Chevallier (2007), we do not inflate all errors by the same factor. For an observation i with reported retrieval error σ i , we consider all N observations within distance R = 500 km and time T = 1 h of i. If all these observations were perfectly correlated, and we averaged them, then the average observation would have an error σ b .
On the other hand, if all N observations were uncorrelated and we inflated their errors by α i and then averaged them, that average would have an error σ uncorr .
If we demand σ b = σ uncorr , then we get the inflation factor α i for the error σ i .
This approach has the benefit that errors of observations that are clustered together are inflated more than of those that are far away from each other, effectively de-weighting clustered observations. This is physically plausible, since errors in observations closer together are expected to be more strongly correlated, resulting in fewer effective constraints. Applied to our X CO 2 dataset, this algorithm results in inflation factors between 1 and 6.

The inversion system
We assimilate surface CO 2 and GOSAT X CO 2 observations in a 4DVAR system to estimate monthly surface fluxes. The inversion system closely follows the framework described in Meirink et al. (2008). We describe it briefly below, and the reader is encouraged to refer to Meirink et al. (2008) for details of the TM5 4DVAR system.

TM5 4DVAR
Let x denote a vector of length N representing surface fluxes, y a vector of length M of atmospheric observations available to constrain the fluxes, and let x 0 be a prior estimate of surface fluxes before any observations are made. We estimate the optimal value of x, given a prior guess x 0 and observations y, by minimizing the cost function.
where H is the composition of an atmospheric transport operator and an observation operator, R is the observation error covariance matrix, and B is the prior flux covariance matrix, whereā denotes the ensemble average of a, i.e., the average over a number of realizations of a. We divide the CO 2 flux into four categories, similar to Peters et al. (2007). The flux from surface grid cell i hor at time t is where the superscripts denote the following: "bio" = flux from terrestrial biosphere, "oce" = flux from oceans, "FF" = fossil fuel emissions, and "fire" = flux from biomass burning and land use change. The flux from each category is specified on a global 6 • × 4 • grid every three hours. x FF and x fire are not optimized; x bio and x oce are optimized on a global 6 • × 4 • grid every month for 22 months. Thus, the size of the state vector is N = 2 × 60 × 45 × 22 = 118 800. Our temporal resolution of one month is coarser than what has been used by other published CO 2 inversions, such as the eight-day resolution of Chevallier et al. (2011) or the weekly resolution of CarbonTracker. However, given our choices of prior temporal correlation described later in Table 1, our inversion system cannot resolve emissions at the sub-monthly scale. Therefore, to avoid unnecessarily inflating the size of the state vector, we solve for fluxes at a monthly temporal resolution.

Prior flux x 0 and flux covariance B
A priori estimates of x bio and x fire are taken from CASA GFED3 (van der Werf et al., 2010). x fire denotes the CO 2 emitted directly due to the burning of vegetation in forest fires, whereas the after effects of the fire -e.g., rotting litter, decaying dead trees -are bundled into the biosphere estimate x bio through heterotrophic respiration (R H ). In the present state of our knowledge, the uncertainty in x fire 0 is much smaller than the uncertainty in x bio 0 , and hence we choose to hold x fire fixed at its prior, inventory-based value, while absorbing small variations in x fire into the much more uncertain x bio . CASA GFED3 provides monthly net primary production (NPP) and R H , from which the net ecosystem exchange, NEE = R H − NPP, is calculated. Higher frequency variations -such as the diurnal cycle -are added to NEE every three hours using a relation involving the two-meter air temperature and incident solar radiation as described in Olsen and Randerson (2004).
A priori estimates of x oce come from an ocean interior inversion described in Jacobson et al. (2007). This inversion yields monthly mean oceanic CO 2 fluxes, which are combined with 3-hourly surface pressures and 10 m winds from the ECMWF ERA 40 model to give synoptic variations and interannual variability as described in Kettle and Merchant (2005).
Fossil fuel emissions x FF , which are held fixed at their prior values in our inversion, come from a combination of several inventories. The procedure for compiling x FF 0 is identical to the fossil fuel module described in the CarbonTracker documentation (CarbonTracker, 2012), except that we do not use the ODIAC database. Annual global total fossil fuel emissions up to 2008 are taken from Boden et al. (2010) and are extrapolated to 2009 and 2010 using BP Statistical Review of World Energy as described by Myhre et al. (2009). The total emission from a country is distributed within the country according to the spatial pattern of the higher resolution EDGAR v4.0 inventory from the European Commission. Where available, a seasonal cycle is imposed on the annual total data, for example by using the monthly data from Blasing et al. (2005). For the two flux categories being optimized, x bio and x oce , we need to specify a spatiotemporal covariance matrix for the prior flux, i.e., the matrix B in Eq. (2). We assume that the covariance is separable in space and time, so that where x r,t is the CO 2 flux from position r at time t. We choose correlation functions C r and C t to be exponentially decaying with characteristic length and timescales L and T , respectively.
Further, we assume no prior correlation between the two optimized categories, biosphere and oceanic fluxes, so each category is characterized by its own L and T . The standard deviation σ r,t is specified as a certain fraction ξ of the absolute flux from (r, t), where σ 0 = 0.005 kg CO 2 /sec/grid box allows the inversion to adjust fluxes over grid boxes with zero prior emission. Each category has its own ξ , and the category-specific parameters are summarized in Table 1. For each category, the parameters in Table 1 were tuned to produce the same ratio of annual flux uncertainty to annual mean flux as seen in the prior fluxes of CarbonTracker 2011, which were 1.20 and 0.51 for the terrestrial and oceanic components, respectively.

Transport and observation operator H
CO 2 surface fluxes were propagated forward in time using the TM5 atmospheric transport model run on a global 6 • ×4 • grid over 25 vertical layers (Krol et al., 2005). We drove TM5 with ECMWF ERA Interim meteorological data, which were coarsened from their original resolution to our model grid. Since the CO 2 mixing ratio in a 3-dimensional grid cell is only representative of the average mixing ratio inside that cell, the mixing ratio at each point sample location was calculated by linear interpolation using the tracer mass and its spatial gradient, both of which are calculated by TM5. A similar linear interpolation in two dimensions was performed when calculating the modeled total column at a satellite sounding location.
For comparing to X CO 2 measurements, the modeled CO 2 column needs to be convolved with the satellite averaging kernel. For this, the column of modeled CO 2 mixing ratio c mod at the time and position of a satellite sounding was first re-binned onto the retrieval grid by multiplying with a redistribution matrix h, then convolved with the averaging kernel A and added to the prior profile c pri to yield the modeled total column where w is a vector of pressure weights, w i = mass of dry air in layer i/mass of dry air in total column, and c pri is the prior CO 2 profile used in the retrieval.

Initial atmospheric CO 2 field
The transport model TM5 simulates atmospheric CO 2 concentrations given a flux scenario and a 3-dimensional CO 2 field at the model starting time, the so-called initial concentration, which is not optimized in our inversions. While small biases due to an incorrect initial concentration can be compensated by surface fluxes during the spin-up period (for us, 1 March-1 June 2009), large biases cannot be corrected in three months and will result in erroneous fluxes over the first few months of our assimilation window. Therefore, it is important to start with an initial concentration field consistent with the state of the atmosphere as captured by our assimilated measurements.
We perform an optimization with only surface CO 2 measurements as described in Sect. 2.1 from 1 January 2008 to 1 May 2009, following the same methodology detailed in this section. The initial concentration for that 4DVAR run is taken from the posterior dry air CO 2 mole fraction field of Carbon-Tracker 2010 on 1 January 2008. The optimized posterior flux over the 16 months is then propagated forward in TM5 to create an optimized atmospheric CO 2 field. This CO 2 field is sampled at 00:00 UTC on 1 March 2009 to create the initial concentration field that is used in our inversions described in the rest of this manuscript. By stopping the 16-month assimilation on 1 May 2009, we make sure that no observation is counted twice, since the inversions presented later assimilate observations made after 1 June 2009.

Observations y and observation error covariance R
Observations used in Eq.
(2) come from one or both of the streams described in Sect. 2. After error inflation (for satellite data) and temporal averaging (for point samples), the observations are assumed to be independent, so R is diagonal, with elements where σ obs,i is the estimated error in observation i (after error inflation, in the case of X CO 2 ), and σ mod,i is the estimated error in modeling the i-th observation. The latter component, which is the so-called model representativeness error, is the error made by our finite resolution model in simulating a sample at a specific location. This error will be large in an area with a large spatial gradient of CO 2 , and will be small if the CO 2 concentration in an area is spatially uniform. For point observations, we estimate this error for a grid cell as the standard deviation in CO 2 mixing ratio across all neighboring grid cells in three dimensions. In case of boundary cells such as surface cells or cells at the poles, the concentration in the non-existent neighbor is simulated by where d is the vector from cell's center to the center of its (missing) neighbor. In the case of point samples, due to the high accuracy of flask and in situ measurements, σ mod,i > σ obs,i for most observations. In the case of satellite measurements, the model representativeness error is calculated for each gridbox (as above) in the modeled column c mod . Let the error in the i-th vertical gridbox in column c mod be e i . Then σ mod for the modeled X CO 2 is calculated by For the TM5 transport model running at 6 • ×4 • lateral resolution and 25 vertical layers, for point samples, the mean σ mod is 3 ppm, larger than the mean σ obs of 2.5 ppm. For satellite observations, the mean σ mod of 0.05 ppm is smaller than the mean σ obs of 3.6 ppm. It should be noted that σ obs = 2.5 ppm for point samples is not much smaller than theσ obs = 3.6 ppm for satellite samples, contrary to the expectation that point samples are much more accurate than satellite samples. This is because the σ obs for satellite samples reflects the estimated error for individual soundings, whereas σ obs for point samples includes the CO 2 variability over a four-hour window that each point sample is supposed to represent. A single point sample is very accurate (estimated error ∼ 0.2 ppm), but over a four-hour averaging window varies by several parts per million, which effectively becomes the margin within which the inversion tries to reproduce them.

The minimization algorithm and the posterior covariance matrix
To minimize the cost function J in Eq.
(2), we use the adjoint TM5 transport model and adjoint point and satellite observation operators to calculate ∇ x J (Kaminski et al., 1999). Details of constructing the adjoint model are given by Meirink et al. (2008) and references therein. Since our problem is linear, i.e., J is a purely quadratic function of x, we use the conjugate gradient algorithm for minimizing J (Navon and Legler, 1987). The application of the algorithm specifically to a 4DVAR data assimilation problem is described by Fisher and Courtier (1995). The method is based on the Lanczos algorithm (Lanczos, 1950) for solving a linear system and allows us to simultaneously minimize J and derive the leading eigenvalues and eigenvectors of the Hessian ∇ 2 x J . Section 2.2 of Meirink et al. (2008) describes how the algorithm is used in a flux inversion system to estimate optimized fluxes as well as the posterior covariance matrix. We summarize it very briefly below.
The optimization is performed in terms of a preconditioned variable z (Courtier et al., 1994).
The Hessian in this preconditioned space is The optimization, after n steps, yields the n leading eigenvalues λ i=1...n and the corresponding eigenvectors v i=1...n of H. These can be used to derive an approximation to the posterior covariance matrixB.
At the start of the optimization, n = 0 and posterior error estimates are as large as prior estimates, orB = B. During an optimization, the λ's start from a high positive value and approach 1 as n → N. Hence with each iteration the posterior covariance matrix is reduced by some amount, but the reduction becomes less and less as λ → 1, untilB reaches its exact value at n = N. In all the inversions reported in this paper, the lowest eigenvalue 1 < λ min < 1.03, and hence our reported flux uncertainties are close to but slightly higher than those estimated from the exact posterior covariance matrix. Uncertainty estimates of spatiotemporally aggregated fluxes, e.g., the global total annual CO 2 emission, can be calculated by summing up the correct elements ofB. Due to the high dimensionality of the problem, exact convergence is never reached, and we stop the conjugate gradient iterations when where η = 10 −8 is a small scalar constant.

Validation of the optimized fluxes
Our flux inversion is an underconstrained problem, i.e., there are fewer measurements than degrees of freedom that need constraining. Therefore, given a prior flux and a set of observations, the optimization will always converge to some optimized posterior flux, which could be quite far from the true flux scenario. A flux inversion should ultimately be judged on how well it can reproduce the atmospheric state of the relevant tracer, which, in our case, is captured by measurements of atmospheric CO 2 concentration. We evaluate our inversions by simulating the atmospheric CO 2 field with posterior flux estimates, and comparing that field with independent CO 2 measurements, i.e., measurements that have not been assimilated by the 4DVAR system. We use two types of measurements, point and total column.

Point samples from CONTRAIL
The Comprehensive Observation Network for TRace gases by AIrLiner (CONTRAIL) project (Machida et al., 2008) has been observing vertical CO 2 profiles over 43 airports worldwide and along intercontinental flight paths using five Japan Airlines (JAL) commercial airliners. The data coverage is extensive in the Northern Hemisphere, and vertically the samples go up to 150 hPa Niwa et al., 2012). We sample our posterior atmospheric CO 2 field at 967 418 points between 1 June 2009 and 1 January 2011 corresponding to CONTRAIL samples available. The sampling algorithm and model error calculation is identical to our point sampling step during assimilation, detailed in Sect. 3.1.4.

X CO 2 samples from TCCON
The Total Carbon Column Observing Network (TCCON) measures X CO 2 at 22 sites globally (Wunch et al., , 2010, a number of which were operational during the time window of our inversion. We sample the posterior atmospheric CO 2 field at 235 467 sampling times and locations between 1 June 2009 and 1 January 2011 corresponding to TCCON samples available. In each case, we convolve our vertical CO 2 profile with the relevant TCCON averaging kernel and add the appropriate TCCON prior CO 2 profile. At most times, the time series of individual X CO 2 soundings has a lot of variability that can obscure the picture if we are interested in phenomena such as the summer drawdown of CO 2 , which happens over several weeks. Therefore, for each site, all the observations and modeled samples have been averaged over seven days to get rid of high frequency variations.

Consistency between surface and GOSAT CO 2
We first check whether there are obvious biases between GOSAT X CO 2 and surface samples by looking at the atmospheric state of CO 2 estimated after assimilating only GOSAT soundings. We compare the posterior CO 2 time series at surface stations in the northern extra-tropics and tropics to observations at those stations in Fig. 3. An inversion using only GOSAT X CO 2 can reproduce the surface time series of CO 2 faithfully, the average model-observation bias de-creasing from a prior value of 1.95 ppm to −0.55 ppm. In the northern extra-tropics, the prior time series do not have sufficiently deep minima in the summer due to too little uptake in the CASA fluxes (Yang et al., 2007). The top row of Fig. 3 shows that assimilating GOSAT X CO 2 corrects for that underestimated uptake. The simulated phasing of the seasonal cycle also matches the observations, even at difficult sites with a very active biosphere such as Park Falls. In the tropics (lower row in Fig. 3), the seasonal cycle is not as prominent, but we still see a GOSAT-only inversion being closer to the observed time series than prior estimates. Although we have shown four stations for brevity, the good agreement between a GOSAT-only inversion and surface stations is prevalent throughout stations in the tropics and the northern extratropics. This positive outcome was by no means obvious at the start of the GOSAT mission, and bodes well for the future of constraining the carbon cycle by remote sensing.
In the southern extra-tropics, however, the picture is not as nice. As shown in Fig. 4, assimilating only GOSAT data overestimates the seasonal cycle. We will see later in Sect. 5 that this stems at least partly from X CO 2 retrieval problems over the ocean, which we are still working on improving. Figure 5 shows monthly averages of CONTRAIL observations in the northern temperate latitudes and the tropics, defined respectively as 23.5-66.5 • N and 23.5 • S-23.5 • N, along with posterior CO 2 concentrations from several inversions sampled at CONTRAIL locations and times. As can be seen in Fig. 5, all inversions produce improved states of the atmosphere compared to the prior. In the northern temperate latitudes, GOSAT data point to a deeper drawdown of CO 2 compared to surface data in the summer of 2009, and a higher outgassing in DJF 2009-2010. We will see in Sect. 4.3 that, as expected, assimilating GOSAT X CO 2 leads to an amplified seasonal cycle in the surface flux. Since the CONTRAIL data do not show a similar amplification -in fact, the CONTRAIL data follow the surface-only inversion quite closely -this enhanced seasonal cycle need not be a real feature of the surface flux. Whether this enhancement is something specific to the 2009-2010 seasonal cycle or is a more general result of assimilating GOSAT data can only be resolved from multi-year GOSAT assimilations in the future.

Comparison to CONTRAIL
In the tropics, assimilating GOSAT data results in a CO 2 drawdown in the Northern Hemisphere (NH) summer of 2009 that is not reflected in the CONTRAIL data. On the other hand, CONTRAIL CO 2 points to an outgassing in the NH spring of 2010 that is better captured by GOSAT than by surface samples. This is also seen in the tropical surface flux time series in Sect. 4.3. It is tempting to investigate whether Atmos. Chem. Phys., 13, 8695-8717 tropics and northern extra-tropics, compared to CO 2 measurements taken at those stations. The "GOSAT" line is from an inversion using only GOSAT X CO 2 , whereas the "Flasks" line is from an inversion using only surface samples. At Park Falls, to prevent the large diurnal cycle in the continuous time series from obscuring the picture, the model has been presented after being co-sampled at observation times and averaged over three-day windows. This is purely for visual clarity. the enhanced seasonal cycle in the southern extra-tropics as seen in Fig. 4 is observed in comparisons to CONTRAIL data as well, but unfortunately the lack of CONTRAIL samples at those latitudes makes such a comparison impossible.

Comparison to TCCON
To evaluate the inversions against TCCON, both TCCON data and co-sampled optimized X CO 2 were averaged over 7day windows to improve visual clarity. The resulting time series are plotted in Fig. 6 (left and center) for the Northern and Southern Hemisphere stations reporting the maximum number of data points over our assimilation time window. Due to data gaps in different TCCON stations over Europe, data from four TCCON stations were combined in Fig. 6 (right). For each station, coincident GOSAT data were chosen based on a sampling procedure by Oshchepkov et al. (2012). First, all GOSAT soundings within ±15 • longitude and ±7.5 • latitude of a station were selected. Of those, soundings were deemed coincident with the TCCON station if the modeled X CO 2 at the sounding locations were within 0.5 ppm of the modeled X CO 2 over the TCCON station. Details of the colocation technique have been discussed by Guerlet et al. (2013). The coincident GOSAT data were averaged over 15 days (owing to higher scatter) and plotted with a star symbol. The shaded region around a GOSAT data point denotes the spread of coincident GOSAT data during an averaging period.
At Lamont, the TCCON data in the summer of 2009 show a drawdown of CO 2 which is mirrored in the GOSAT data and in inversions with GOSAT data, but not in inversions with surface data (Fig. 6, left). This points to enhancements of the summer uptake over the continental United States or upwind areas that are not captured by surface CO 2 data.  Fig. 4. Model estimated CO 2 time series at two surface stations in the southern extra-tropics, compared to CO 2 measurements taken at those stations. The "GOSAT" line is from an inversion using only GOSAT X CO 2 , whereas the "Flasks" line is from an inversion using only surface samples.
In the future, with the availability of multi-year GOSAT datasets, we will investigate whether 2009 had anomalously high uptake or whether the mismatch between surface and GOSAT inversions persists year after year.
A similar enhancement is not seen at European TCCON stations; an inversion of only surface data reproduces the summer 2009 drawdown seen by TCCON stations quite well (Fig. 6, right). Thus the enhanced depletion of summer 2009 CO 2 higher up, which is sensed by GOSAT over North American TCCON stations, is a regional/continental-scale phenomenon rather than one occurring across the entire latitude band. This is consistent with the fact that CONTRAIL samples over the northern temperate latitudes (Fig. 5, left) do not see an enhanced depletion of CO 2 when compared with a surface-only inversion.
At Wollongong, the bias between TCCON and co-located GOSAT data is seen to vary seasonally. Apart from a small overall high bias in coincident GOSAT X CO 2 , the GOSAT data in July-August 2010 is higher than TCCON by ∼ 2 ppm, whereas in March-April 2009 it is higher by ∼ 1 ppm. This seasonal variation in the bias enhances the seasonal cycle of the inverted fluxes in the southern temperate latitudes, which is consistent with the enhanced seasonal cycle of CO 2 seen at surface stations in Fig. 4.

Annually integrated fluxes
In Figs. 7 and 8, the optimized fluxes are presented after aggregation over one year and different land and ocean regions. For reference, the regions are defined as global 0.5 • × 0.5 • masks in the Supplement. The error bars are 1σ posterior errors. The rightmost group of bars in Fig. 8 denote the global budget of CO 2 . The difference between the global budgets predicted by a surface-only inversion and a joint surface + GOSAT inversion corresponds to a difference of 0.5 ppm in the global average CO 2 mixing ratio. From an 18-month inversion, it is not possible to constrain the global growth rate to a higher accuracy, and future inversions over longer periods will be able to separate the trend from interannual variability. Of the additional 1 PgC predicted by the joint inversion, 90 % lies above 850 hPa, which is consistent with the fact that in the surface layer both the surface-only inversion and the joint inversion are consistent with surface stations. Figure 9 shows the posterior correlation matrix between the different emissions shown in Figs. 7 and 8, for the inversion assimilating both surface CO 2 and GOSAT X CO 2 measurements. For geographically disjoint regions, a near-zero correlation would imply that the emissions from those regions can be independently constrained by the observations, whereas a high negative correlation would indicate a compensatory mechanism. To conclude that optimized emissions over a certain region during a certain period are robustas opposed to artifacts of inversion -there must be significant uncertainty reduction, i.e., in Figs. 7 and 8 the posterior error bar must be smaller than the prior error bar. Moreover, the correction to prior emission over that region must not be strongly negatively correlated to the correction over a geographically disjoint region in Fig. 9. This would mean, for example, that although the global 2.3 PgC sink between 1 September 2009 and 1 September 2010 is a robust conclusion from a joint GOSAT + surface CO 2 inversion, its partitioning between a 2 PgC land source and a 4.3 PgC ocean sink is not; a correlation of −0.97 between global land and global ocean in Fig. 9 points to an inability of the inversion to independently constrain the emissions from those two regions. Similarly, although tropical Asia does not have a significant correlation with any other TRANSCOM region in Fig. 9, it shows only 15 % uncertainty reduction in a joint inversion, limiting the trustworthiness of the estimated 0.2 PgC source.

Seasonal variation of CO 2 flux
We group the optimized flux over 18 months for each region into blocks of three months, which splits a year into summer, fall, winter and spring. Figure 10 (top row) shows that in the northern extra-tropics, assimilation of GOSAT data leads to a slightly higher uptake in the NH summer of 2009 (JJA '09) and slightly higher outgassing in the NH winter of 2009-2010 (DJF '09). The strengthening of uptake and outgassing are independent events and not compensations for each other to maintain a consistent annual budget, since they have no posterior correlation (Fig. 11, left). Since CONTRAIL data in Sect. 4.2.1 do not support this enhancement, it could be an artifact of transport error as opposed to a real signal. Whether this is specific to the seasonal cycle in 2009 or whether this is a more generic feature can only be answered in the future with multi-year GOSAT inversions. Figure 10 (middle row) shows that in the tropics assimilating GOSAT X CO 2 results in a CO 2 drawdown in the NH summer of 2009 due to increased uptake by the tropical land, and increased outgassing from the tropical land in the NH spring of 2010 leads to a CO 2 increase. The two events are independent (Fig. 11, center), and it is likely that the increased outgassing in MAM 2010 is a genuine surface flux signal, since it is also visible in CONTRAIL samples over the tropics (Fig. 5, right).
Assimilating GOSAT X CO 2 increases the seasonal cycle amplitude of the CO 2 flux from the southern extra-tropics, as seen in Fig. 10 (bottom row). In Southern Hemisphere (SH) spring (SON '09 and SON '10) GOSAT X CO 2 points to an outgassing of CO 2 from the land, while in SH summer (DJF '09) GOSAT X CO 2 indicates a deeper oceanic sink of CO 2 than what is consistent with surface measurements. This leads to a net amplification of the CO 2 seasonal cycle. This enhancement is seen neither by the southern temperate TC-CON stations such as Wollongong (Fig. 6), nor by surface stations in those latitudes such as Crozet Islands (Fig. 3). Since this enhanced seasonality is present in the GOSAT X CO 2 as evidenced by its comparison with TCCON X CO 2 at Wollongong in Fig. 6 (center), there is likely a retrieval artifact over the southern temperate latitudes that leads to a spurious strengthening of the CO 2 seasonal cycle in the south. Alternatively, there could also be a retrieval artifact over oceans, which would alias into the X CO 2 seasonal cycle over the southern temperate latitudes owing to the high percentage of ocean cover there and their seasonally varying satellite coverage. Indeed, we will see in Sect. 4.5 that allowing for a global land-sea bias in the assimilated X CO 2 and estimating it during the inversion gets rid of part of the spurious seasonal cycle enhancement seen in Fig. 4.

Robustness of optimized fluxes
There are many sources of error in the 4DVAR inversion system, of which only two -namely, error in the specification of prior fluxes (B in Eq. 2) and error in simulating observations (R in Eq. 2) -are explicitly incorporated in the data assimilation system. Our final flux estimates and the associated error bounds reflect the uncertainties introduced by these two error sources. There are other, more systematic errors, however, that are not accounted for byB in Eq. (10). The most critical of them are systematic errors in atmospheric tracer transport, of which there could be many sources, such as finite spatiotemporal resolution, incorrect vertical transport, systematic errors in the driving meteorological fields, and incorrect interhemispheric exchange rates. Estimated fluxes are also sensitive to the form of spatiotemporal covariance assumed  Fig. 6. 7-day averaged TCCON data at Lamont (left) and Wollongong (center), along with posterior CO 2 fields co-sampled and convolved with the respective averaging kernels. For each TCCON station, coincident GOSAT observations have been averaged over 15 days -due to higher scatter -and plotted over. The shaded area is the standard deviation of all the coincident data for a month, i.e., it is the "error bar" region of the stars. Over Europe (right), due to gaps in TCCON data over individual stations, data from four stations have been combined, along with their corresponding coincident GOSAT data and posterior co-sampled X CO 2 .   Table 1. Section S1 of the Supplement describes our sensitivity tests exploring the robustness of estimated fluxes to the factors mentioned before. We see that our estimated fluxes, aggregated annually over continental length scales, are robust, since none of the tests we performed could change the aggregates beyond 1σ .

Optimizing a land-ocean bias
Figures 4 and 6 show that assimilating GOSAT X CO 2 leads to an overestimation of the seasonal cycle in the southern extra-tropics. We also know that the land-ocean partitioning of the global carbon sink in Fig. 8 is not realistic, because it shows a large net land source of carbon, contrary to the present understanding of the terrestrial carbon cycle. We suspect that the two problems have a common origin, viz. problematic X CO 2 retrievals over oceans. Since most of the southern extra-tropics is covered by oceans, which has a seasonally varying coverage from GOSAT due to the strict requirements of the glint geometry, a systematic land-sea bias in GOSAT X CO 2 would alias into the seasonal cycle of the southern extra-tropics. To investigate the plausibility of a land-sea bias in our GOSAT X CO 2 retrievals, we introduce two bias parameters, b land and b ocean , with the interpretation that retrieved values of X CO 2 are off from the "true" values of X CO 2 by b land and b ocean over land and ocean, respectively. X CO 2 (land,true) = X CO 2 (land,retrieved) + b land X CO 2 (ocean,true) = X CO 2 (ocean,retrieved) + b ocean (12) We then try to find the optimal values of b land,ocean by an inversion assimilating GOSAT and surface layer CO 2 data. We start with prior values of 0 ppm for both, and since land pixels are better calibrated than ocean pixels owing to TCCON stations on land, we assign prior errors of σ (b land ) = 1 ppm and σ (b ocean ) = 2 ppm. Even though we are primarily interested in a land-ocean bias, we optimize not one but two  parameters, b land and b ocean , since a single degree of freedom would imply that the overall level is perfectly constrained, which is not the case. At the same time, the offsets over land and ocean should not be completely independent, so we constrain b to have 1.4 degrees of freedom by assigning a prior correlation of corr(b land , b ocean ) = −0.65. After assimilating GOSAT and surface CO 2 data, the optimal values of the bias parameters are b land = −0.14 ± 0.04 ppm and b ocean = 0.79 ± 0.04 ppm. As seen in Figs. 12 and 13, with this simple bias correction, the optimized fluxes (labeled "(BC)" in subsequent figures) have much more realistic land-ocean partitioning. The global land is now a 0.6 ± 0.4 PgC yr −1 sink instead of a net carbon source, and the global ocean is now a 1.7 ± 0.4 PgC yr −1 sink, which is very similar to a surface-only inversion. The tropical ocean becomes a 0.8 ± 0.2 PgC yr −1 source, and the poleward carbon In the northern extra-tropics, since the seasonal cycle imposed by the terrestrial biosphere is much larger than that imposed by the ocean, almost all of the northern extra-tropical seasonal signal comes from the land. In the tropics and southern extra-tropics, the magnitude of the seasonal variation of the oceanic source/sink is comparable to that of the land source/sink, although the phasing is different. Atmos. Chem. Phys., 13, 8695-8717, 2013 www.atmos-chem-phys.net/13/8695/2013/ flux, i.e., the carbon sourcing in the tropics and sinking in the extra-tropics, is increased compared to a surface-only inversion.
Optimizing for a land-ocean bias in GOSAT X CO 2 in our data assimilation also reduces the enhancement of the seasonal cycle in the southern extra-tropics discussed before. Figure 14 shows that if we account for a possible land-ocean bias in assimilated X CO 2 , the posterior CO 2 time series at the surface in the southern extra-tropics displays a higher concentration of CO 2 in the NH winter of 2009-2010. This reduces the spurious enhancement of the seasonal cycle at these latitudes as mentioned before. We conclude that small biases in our ocean retrievals, coupled with the seasonally varying coverage of GOSAT ocean pixels, results in the spurious enhancement of the southern extra-tropical CO 2 seasonal cycle. Ocean pixels have a much higher signal to noise ratio due to the high surface reflectivity in sunglint geometry, so in principle they could be more accurate than land pixels. Therefore, we are currently working on improving our ocean X CO 2 retrievals to eliminate the remaining biases. Pending further improvements to our ocean retrievals, we consider a joint inversion of surface and satellite data along with the optimization of two bias parameters from Eq. (12) as the most reliable method for assimilating RemoTeC X CO 2 in a sourcesink inversion. Below, we briefly summarize optimized emissions from this joint inversion where the numbers quoted hold up to scrutiny both in terms of uncertainty reduction and posterior correlation, as mentioned in Sect. 4.3.1.
A surface-only inversion changes the 0.4 ± 0.5 PgC North American (prior) source into a 0.4 ± 0.2 PgC sink. Adding GOSAT X CO 2 strengthens this sink to 1.0 ± 0.1 PgC. Europe goes from a 0.3 ± 0.4 PgC prior source to a 0.3 ± 0.3 PgC sink in a surface-only inversion, with the addition of GOSAT strengthening the sink to 1.3 ± 0.2 PgC, although part of it is to compensate the sourcing from the North Atlantic. South America, a 0.5 ± 0.5 PgC sink in the prior, is slightly weakened by the surface inversion to 0.4 ± 0.4 PgC, and assimilating GOSAT X CO 2 turns it almost carbon neutral, a small 0.1±0.2 PgC source. The Eurasian boreal region, which was a 1.1 ± 0.3 PgC sink according to a surface-only inversion, is considerably weakened to 0.3 ± 0.2 PgC on assimilating GOSAT. The Eurasian temperate region is turned from a 0.1 ± 0.2 PgC sink (surface-only) to a 0.3 ± 0.2 PgC source (surface + GOSAT), turning Asia into a net source of 0.3 ± 0.3 PgC from a surface-only estimate of 1 ± 0.4 PgC sink. Part of this, however, could be a compensatory effect to offset the enhanced sinking in the North Pacific. Overall, the prior 0.3 ± 0.7 PgC sourcing from the tropics is strengthened by a surface inversion to 0.5 ± 0.4 PgC, which is further enhanced to a 2.1 ± 0.2 PgC source on addition of GOSAT X CO 2 . The land-ocean partitioning of this source (1.3 ± 0.2 PgC from land and 0.8 ± 0.2 from ocean) is less trustworthy owing to a −0.57 posterior correlation between tropical land and tropical ocean. A part of this additional CO 2 is taken up by the northern extra-tropical sink (2.1 ± 0.3 PgC for a surface-only inversion compared to 2.9 ± 0.1 PgC for a surface + GOSAT inversion), and the rest contributes to the weakening of the global sink from 3.4±0.2 PgC (surface only) to 2.3±0.1 PgC (surface + GOSAT).

Discussion
In Sect. 4 we presented the results of our source-sink inversions between 1 June 2009 and 1 December 2010 using both GOSAT X CO 2 and surface measurements of CO 2 . We saw in Sect. 4.1 that when we optimized surface CO 2 fluxes against only GOSAT X CO 2 , we could reproduce the CO 2 time series at surface stations in the tropics and northern extra-tropics. In the northern extra-tropics, assimilating only GOSAT X CO 2 corrects the underestimated summer uptake of the CASA prior fluxes, and also improves the phasing of the seasonal cycle. Surface fluxes optimized against GOSAT measurements could also reproduce the X CO 2 time series at TCCON stations situated in these latitudes, evidenced in Fig. 6 at Lamont and several European stations. We stress that we have demonstrated that X CO 2 measurements from GOSAT and similar future satellite missions can impose significant quantitative constraints on the carbon cycle, a result that bodes well for future remote sensing missions. However, our results also point to remaining problems in Re-moTeC retrievals over the southern extra-tropics. Inversions assimilating GOSAT X CO 2 show enhanced seasonal cycles not seen in other independent atmospheric measurements. Corresponding inversion-derived seasonal emission adjustments are found over the southern extra-tropics, such as CO 2 emissions from the southern extra-tropical land during SON '09 and SON '10, and enhanced CO 2 uptake by the southern extra-tropical ocean in DJF '09 ( Fig. 10, bottom row). The effect of this enhancement is visible in the simulated surface layer CO 2 time series at surface stations in the southern extra-tropical latitudes, such as Crozet Islands in Fig. 3. This enhancement is also visible in simulated TCCON time series from those latitudes, such as at Wollongong in Fig. 6. Since neither the surface stations nor TCCON see this enhanced seasonal cycle, we suspect that this enhancement is a spurious effect of retrieval artifacts as opposed to a true surface flux signal. In fact, co-located GOSAT data over Wollongong (Fig. 6, center), compared to TCCON X CO 2 , show a higher bias in SON '09 and SON '10 as opposed to MAM '10, indicating that this enhancement could be a direct effect of the X CO 2 data ingested instead of a transport model error. To be certain of this conclusion, we need to compare GOSAT X CO 2 in this latitude band with ground-based X CO 2 measurements from several locations. Unfortunately, the only other TCCON station in the southern extra-tropics, Lauder, has a data gap during MAM '10.
Following the method of Conway et al. (1994), the increase in the global average CO 2 mixing ratio derived by NOAA between 1 September 2009 and 1 September 2010    is 2.72 ppm (data downloaded from ftp://ftp.cmdl.noaa.gov/ ccg/co2/trends/co2 mm gl.txt). The global CO 2 emission over that period predicted by our inversions (Fig. 8) corresponds to increases of 2.4 ± 0.1 ppm (surface data only), 2.99 ± 0.04 ppm (GOSAT only) and 2.86 ± 0.04 ppm (joint assimilation). All three estimates are in line with the NOAA estimate, and the remaining discrepancies can be explained by the fact that while the NOAA estimate is the result of a multi-year trend analysis, we have only optimized the flux for 18 months, and are thus significantly affected by interannual variability. The trend estimation method of Thoning et al. (1989) used by Conway et al. (1994) cannot be used in our case since the data record is too short. In the future, with multi-year GOSAT inversions, the comparison with the NOAA flask-derived estimate of CO 2 growth rate could prove to be more insightful.
Although the global CO 2 budget in Fig. 8 is consistent with independent estimates such as the one from NOAA, the land-ocean partitioning of the global total flux is less reliable; in particular the global land biosphere turning out to be a net source of carbon is not believable. We suspect that this erroneous partitioning is due to a retrieval problem over the ocean, since land retrievals are rather well calibrated against TCCON measurements. We show that optimizing for a global bias between land and ocean X CO 2 retrievals not only makes the partitioning more realistic, it also reduces the spurious enhancement of the seasonal cycle seen in the southern extra-tropics. The performance of a joint inversion which optimizes a global land-sea bias in X CO 2 is good enough for us to consider it as our standard joint inversion.
We would like to stress here that a total land-ocean bias of 0.93 ppm was enough to induce the dramatic changes in  Fig. 4. The inversion labeled "+ GOSAT (BC)" is a joint inversion in which a land-ocean bias in GOSAT X CO 2 data has been optimized for, whereas the inversion labeled "+ GOSAT" is a joint inversion without any bias correction. The optimized CO 2 time series at the surface in the SH extra-tropics from "+ GOSAT (BC)" has a lower seasonal cycle amplitude compared to "+ GOSAT".
land-ocean partitioning of the fluxes seen in Fig. 13. This emphasizes the stringent accuracy requirements on satellite X CO 2 to be used for source-sink inversions of CO 2 . A 0.93 ppm bias between land and ocean pixels is not something that can be constrained by the present GOSAT X CO 2 validation technique of comparison to X CO 2 from the TC-CON network, but such a bias is enough to significantly change the picture of the global carbon cycle arrived at through an atmospheric inversion. Though the TCCON network is in theory accurate enough to detect this 0.93 ppm bias, the current lack of marine TCCON instruments makes this detection impossible. We would also like to caution that studies assimilating any total column data in the near infrared -be it TCCON, GOSAT or some other satellite -could in principle suffer from fair weather bias, since only X CO 2 retrievals from scenes taken during clear sky and low aerosol conditions are assimilated. The impact of this would be felt in the tropics due to intermittent cloud cover, and in the northern high latitudes in NH winter for similar reasons. Areas of high aerosol content, such as downwind of the Sahara and the Gobi deserts, could also be undersampled, leading to biases.
Based on surface CO 2 time series and source-sink inversions assimilating surface CO 2 observations, previous studies have estimated the northern extra-tropics to be a strong sink of carbon, whereas the tropics have been estimated to be a strong source (Gurney et al., 2002;Gurney, 2004;Tans et al., 1990). The magnitude of this carbon flux, however, has been the subject of some debate. Stephens et al. (2007) showed that global flux scenarios that had weaker poleward carbon fluxes produced atmospheric CO 2 fields more consistent with aircraft measurements of CO 2 . We could thus expect that incorporating CO 2 measurements from aircrafts in a source-sink analysis would weaken the strengths of the tropical source and extra-tropical sink, compared to inversions assimilating only surface CO 2 data. In this study, we find that assimilating GOSAT X CO 2 in addition to surface CO 2 data strengthens the poleward carbon flux, i.e., the tropics become a larger source and the extra-tropics become a deeper sink. As shown in Fig. 13, in a surface data-only inversion, the tropics are a net source of 0.5 ± 0.4 PgC yr −1 , whereas the northern extra-tropics are a net sink of 2 ± 0.3 PgC yr −1 . After assimilating GOSAT X CO 2 , the tropics source 2.1±0.2 PgC yr −1 , and the northern extra-tropics sink 2.9±0.1 PgC yr −1 . We believe this enhancement of the poleward carbon flux to be a robust result within our inversion framework, since it survives all our sensitivity tests performed in Sect. S1 of the Supplement. We note, however, that Stephens et al. (2007) found a large spread in the estimate of the northern extra-tropical sink between different inversions due to differences in vertical transport, and we need to assess the accuracy of vertical transport in TM5 by comparing modeled CO 2 with aircraft profiles to have more confidence in our result. Further, this enhancement of the tropical source and extra-tropical sink, even if robust, could be a feature of the period presented in Sect. 4.3.1, and in the future we plan to explore this question further with multi-year inversions.
The majority of the northern extra-tropical carbon sinking mentioned by Tans et al. (1990) and others happens in the forests and grasslands of North America and boreal Eurasia. With globally increasing CO 2 levels and the changing climate, it is an open question whether these sinks are becoming more or less powerful. Existing biosphere models underestimate the size of these sinks (Messerschmidt et Fig. 15. Comparison of our optimized fluxes with those from CarbonTracker 2011 over the TRANSCOM regions and the global land and ocean. All fluxes are without fossil fuel emissions. Since TRANSCOM regions include northern and southern Africa, whereas we aggregate over Saharan and sub-Saharan Africa (Fig. 7), in this figure we have presented the flux total over Africa.
2012; Yang et al., 2007), as do some inversions using only surface CO 2 measurements (Basu et al., 2011). Therefore it is important to impose additional observational constraints to estimate these sinks more accurately, and then to monitor the evolution of these sinks over time. Our estimate of the North American carbon sink (1.0±0.1 PgC yr −1 ) is higher when we assimilate GOSAT X CO 2 and surface CO 2 as opposed to when we assimilate only surface CO 2 (0.4 ± 0.2 PgC yr −1 ). We believe this result to be robust within our inversion framework and not an artifact of inversion since Fig. 9 shows no appreciable posterior correlation between North America and other regions. For the same reason, the weakening of the Eurasian boreal sink (−1.1 ± 0.3 PgC yr −1 for a surface-only inversion compared to 0.3±0.2 PgC yr −1 when GOSAT data are included) could be a robust result as well. An enhanced uptake over North America and a reduced uptake over Eurasia fits the observed scenario that the TCCON time series at North American sites such as Lamont (Fig. 6) show a summer uptake beyond that predicted by surface data, whereas those over European sites do not. CONTRAIL samples over this latitude band, which are more indicative of the uptake over the entire northern temperate latitudes rather than individual continents, also do not show an enhanced uptake (Fig. 5). Whether these are specific to the summer of 2010 or systematic remains to be seen from future multi-year inversions with GOSAT data. It is interesting to compare the flux estimates we obtain from a joint inversion with bias optimization for the period 1 September 2009 to 1 September 2010 with other inversion products. There are only a few top-down estimates for this period, and we choose two for comparison, (a) the recently released CarbonTracker 2011 (CT2011) and (b) an atmospheric inversion assimilating TCCON X CO 2 by Chevallier et al. (2011). Figure 15 compares CT2011 fluxes with aggregated fluxes from our surface data-only and joint inversions for the TRANSCOM regions and global land and sea. We assimilate the same surface dataset as CT2011 and use the same tracer transport model (TM5), and therefore we could expect our surface-data only inversion to correspond closely to CT2011. It turns out that aggregated fluxes from our surfaceonly inversion and CT2011 fluxes are within 1σ of each other, except over parts of the Pacific. Comparing CT2011 fluxes to our standard product, i.e., a joint inversion with optimized land-sea bias correction, turns up more significant differences. At the large scale, the most significant difference between CT2011 and our bias-optimized joint inversion is that our inversion predicts a significantly smaller land sink (0.6 ± 0.4 PgC) compared to CT2011 (1.9 PgC). The difference does not stem from a single region, but is spread over multiple areas. The North American and European sinks are stronger in our inversion, but the Eurasian boreal sink is weaker. The Eurasian temperate sink seen in CT2011 is a source for us, and the African source of CT2011 is enhanced significantly. We caution that GOSAT X CO 2 over large parts of Africa and the Eurasian temperate region cannot be calibrated against any TCCON measurements, while the Eurasian boreal region has seasonal coverage and is therefore subject to a seasonal sampling bias. Hence the flux estimates over these three regions from our joint inversion should be interpreted with caution.
The global oceanic sink is the same for both CT2011 and our joint inversion, but the distribution is not identical. The Southern Ocean is a stronger sink in our joint inversion, as are the North Pacific temperate and the South Atlantic temperate regions. On the other hand, the South Pacific temperate region is a weaker sink. The North Atlantic temperate region is a source in our inversion, which does not agree with our existing knowledge of the oceanic carbon cycle, and  therefore points to remaining retrieval issues not addressed by our simple bias correction scheme. According to Chevallier et al. (2011), TCCON X CO 2 mostly constrains terrestrial sources and sinks, so in Chevallier et al. (2011) they presented aggregated annual fluxes over the 11 TRANSCOM land regions. We compare our estimates between 1 September 2009 and 1 September 2010 with their estimates over the same period after deducting fossil fuel fluxes in Fig. 16 (F. Chevallier, personal communication, 2012). The global land sink from our inversion is weaker compared to Chevallier et al. (2011), just as it was weaker compared to CT2011. The differences over land are most striking in the Americas. We predict a much stronger carbon sink over North America and a small carbon source over South America, compared to a net South American sink in Chevallier et al. (2011). Interestingly, the posterior atmospheric CO 2 fields from both inversions are consistent with TCCON X CO 2 measurements over North America, because Chevallier et al. (2011) assimilated TCCON X CO 2 while we assimilated GOSAT X CO 2 which were validated against TC-CON X CO 2 . This suggests that at present, different X CO 2 measurements consistent with the same set of TCCON X CO 2 can yield dramatically different posterior flux distributions, even over regions such as North America, which has several TCCON stations.

Conclusions
In this manuscript we have presented optimized global source-sink estimates of CO 2 using the RemoTeC retrieval of GOSAT X CO 2 over eighteen months from 1 June 2009 to 1 December 2010. We have compared the flux estimates with a more conventional 4DVAR inversion using data only from surface-layer CO 2 measurements. We have shown that GOSAT X CO 2 is consistent with surface stations in the tropics and northern extra-tropics, whereas in the southern extratropics it has a higher seasonal cycle compared to surface stations. We believe that this enhanced seasonal cycle is due to a retrieval artifact over the oceans. We show that optimizing a global bias between land and ocean retrievals of GOSAT X CO 2 points to a 0.93 ppm remaining bias between them, and demonstrate that this sub-ppm bias significantly changes the picture of the terrestrial carbon cycle we can derive by assimilating GOSAT X CO 2 . It is worth noting in this context that GOSAT is a new instrument, and as we -as a community -keep working with GOSAT data, we will learn about such biases, their effects, and how to correct them so as to use GOSAT to make concrete statements about the terrestrial carbon cycle. Therefore, our work presented in this manuscript should be considered not just as a collection of new insights into the carbon cycle, but more importantly as a demonstration of the capability of GOSAT to constrain the carbon cycle. Those capabilities should improve when we apply the method described here to future CO 2 -sensing missions such as OCO-2 and CarbonSat.
We have shown that the global budget is well-constrained by GOSAT X CO 2 to yield a global CO 2 growth rate consistent with the NOAA estimate using the method of Conway et al. (1994), although the land-sea partitioning of the global flux remains problematic. We find that both the tropical source and the extra-tropical sink of carbon estimated in earlier works (Gurney et al., 2002;Gurney, 2004;Tans et al., 1990) are strengthened when we include GOSAT X CO 2 in our inversion. This is a robust result within our inversion framework, although it contrasts with earlier findings by Stephens et al. (2007). We also find that the North American carbon sink is slightly stronger than predicted by a surfacedata-only inversion, although North America is well covered by the surface measurement network and the effect of introducing GOSAT data is small.
We caution that these are results using one retrieval (Re-moTeC) for slightly more than one year, and the RemoTeC L2 data product is under constant development. More importantly, our estimate of transport uncertainty is likely to be too small since we have used a single transport model. Although we did perform several sensitivity tests as shown in Sect. S1 of the Supplement, it remains to be seen whether our conclusions hold up when a different -possibly improved -retrieval dataset is assimilated or when a different transport model is used, and if the features we observe are specific to our inversion period or if they recur every year. As the GOSAT X CO 2 dataset is extended through 2011 and 2012, we hope to address those questions in our future work with multiple GOSAT X CO 2 retrievals and over multiple years.