Using continental observations in global atmospheric inversions of CO2: North American carbon sources and sinks

We evaluate North American carbon fluxes using a monthly global Bayesian synthesis inversion that includes wellcalibrated carbon dioxide concentrations measured at continental flux towers. We employ the NASA Parametrized Chemistry Tracer Model (PCTM) for atmospheric transport and a TransCom-style inversion with subcontinental resolution. We subsample carbon dioxide time series at four North American flux tower sites for mid-day hours to ensure sampling of a deep, well-mixed atmospheric boundary layer. The addition of these flux tower sites to a global network reduces North America mean annual flux uncertainty for 2001–2003 by 20% to 0.4 Pg C yr-1 compared to a network without the tower sites. North American flux is estimated to be a net sink of 1.2 ± 0.4 Pg C yr-1 which is within the uncertainty bounds of the result without the towers. Uncertainty reduction is found to be local to the regions within North America where the flux towers are located, and including the towers reduces covariances between regions within NorthAmerica. Mid-day carbon dioxide observations from flux towers provide a viable means of increasing continental observation density and reducing the uncertainty of regional carbon flux estimates in atmospheric inversions.


Introduction
About half of the anthropogenic carbon emitted into the atmosphere remains in the atmosphere each year. The remainder is taken up by the ocean and terrestrial ecosystems through the processes responsible for the natural exchange of carbon between the atmosphere and terrestrial vegetation and the surface ocean (Denman et al., 2007;Forster et al., 2007). Numerous studies (Myneni et al., 2001;Nemani et al., 2003;Potter et al., 2003) show that climate cycles, local weather and ecosystem conditions all affect the interannual variability of this uptake of carbon. Our understanding of the mechanisms governing the dynamics of the carbon cycle has been hampered by a limited ability to locate and quantify these exchanges at sufficiently fine temporal and spatial resolution (Bousquet et al., 2000;Gurney et al., 2002;Ciais et al., 2005;Baker et al., 2006;Peters et al., 2007). Accurate and precise quantification of sources and sinks at regional and continental scales is likely to be increasingly important for evaluation and monitoring of carbon management policies.
Global atmospheric inversions have been used to infer sources and sinks of carbon (both natural and anthropogenic) at continental and ocean basin scale from atmospheric measurements of carbon dioxide using tracer transport models. Model intercomparison projects, including the TransCom (Atmospheric Tracer Transport Model Intercomparison Project) series Gurney et al., 2004;Baker et al. 2006), have been designed to attribute the uncertainties in the continental and ocean basin fluxes estimated by this method. These studies show that transport model differences and the uneven and sparse global distribution of atmospheric carbon dioxide measurements contribute to the uncertainty of the inverse flux estimates. While transport models are improving and the global measurement network for carbon dioxide is expanding, there are still fundamental representation and aggregation errors (Kaminski et al., 2001;Engelen et al., 2002) inherent in the global atmospheric inversion method. There is a mismatch in space and time resolution between the transport models (grid boxes and minutes), the observations (points in space and time) and the inversion solution (continents or subcontinents and months or weeks). Observations are subject to local atmospheric variations. These subgrid scale mesoscale variations are not explicitly accounted for in the transport models; however, we use these local observations to constrain continental and ocean basin results. In addition, global inversions typically require strong and uncertain assumptions about the correlation of fluxes and observations in space and time. If the inversion solution is constructed at the continental scale, for example, it is not possible to evaluate changes in fluxes from subregions within the continent. These strong assumptions about correlations of fluxes in space and time are a weakness; the strengths of such a continental scale global inversion method are the minimum number of unknowns and the global coverage. Some assumptions about coherence in space and time are essential; atmospheric observations will always be uneven and sparse at some level of resolution.
A compelling approach is to invert on the grid and time scale of the transport model (Kaminski et al., 2001;Engelen et al., 2002). Global atmospheric inversions at the spatial resolution of the transport model (e.g. Rödenbeck et al., 2003a;Gourdji et al., 2008;Mueller et al., 2008), aim for the finest resolution possible to minimize representation assumptions, at the expense of larger posterior covariances. Subsequent aggregation into coarser regional and temporal resolution is then used to lessen the posterior error. Regional atmospheric inversions target a geographically limited domain with finer spatial and temporal resolution (Gerbig et al., 2003;Peylin et al., 2005;Lauvaux et al., 2008;Schuh et al., 2009). Both of these approaches involve many more unknowns, which cannot be resolved independently given the current observation density. The underlying assumptions may be minimized, but at the expense of building prior covariance matrices and the increased computational costs required by the finer resolution.
In this experiment we take a pragmatic, middle-ground approach to the continental-scale global inversion by choosing a number of regions roughly matched to the observation density currently available. If we have chosen observation sites that are representative of the regions and sensitive to the surface exchanges in these regions, then we expect that posterior uncertainties and spatial correlations will be reduced and that the problem will be computationally tractable using simple inversion methods. Inversion results can be aggregated to the larger TransCom continental regions for comparison with published results. We can also test the ability of the expanded network to constrain the smaller regions with this method.
Typically global atmospheric measurement network sites have been chosen to facilitate sampling background concentrations of trace gases including carbon dioxide. These background measurement networks have yielded important understanding of interhemispheric gradients in carbon dioxide mixing ratios (Tans et al., 1990;Denning et al., 1995;Keeling et al., 1996) and of the mean annual cycles of carbon emissions and uptake (e.g. Keeling et al., 1995). These data, however, provide limited understanding of the continental carbon cycle. We cannot diagnose continental or regional scale fluxes and determine the factors influencing terrestrial fluxes without observing sites over the continents. Continental carbon dioxide measurements are characterized by strong diurnal and seasonal cycles that reflect a combination of biological fluxes and atmospheric boundary layer dynamics (Bakwin et al., 1998;Yi et al., 2001;Davis et al., 2003). Continental data also contain strong gradients driven by weather (e.g. Hurwitz et al., 2004;Wang et al., 2007;Parazoo et al., 2008). These strong, rapidly varying gradients in the observations may be difficult to simulate in the transport, but the continental data contain information needed to resolve regional sources and sinks of carbon with increasing spatial and temporal resolution.
In this paper, we use the Bayesian synthesis inversion method to demonstrate the impact of including more continental measurement sites in the global measurement network. The added sites are long-running eddy covariance flux towers with high precision carbon dioxide measurements calibrated to global standards. Carbon dioxide measurements at flux towers do not need to be calibrated to global standards for the calculation of net ecosystem exchange of carbon dioxide using the eddy covariance method. The sites used in this study, however, are part of a growing network where the calibration is done with the intent of providing data suitable for application to atmospheric inversion studies. The five towers used in this study have data available during the 2000-2004 time period. Increasing numbers of flux towers are incorporating the calibration processes into their routine processing; this offers opportunities for extending this research in the future.
We focus here on the effect on the North American carbon balance, recognizing the danger that, in an ill-conditioned problem such as this, increasing the density of observations in North America may introduce new challenges. For example, global inversions typically exhibit dipole behaviour or 'pair-sum' relationships (Rödenbeck et al., 2003a) where the flux as a sum for two regions can be constrained, while the individual regions cannot. Within North America we may discover dipoles between the subregions that were not apparent in the continentalscale inversion. With the exception of Boreal Asia, the Northern Hemisphere is well represented in the networks tested (Fig. 2). Concentrating observation sites in North America may highlight dipole relationships between Boreal Asia and other regions of the Northern Hemisphere. Some recent inversions, for example, find a larger terrestrial carbon sink in Europe (Baker et al., 2006;Mueller et al., 2008) and others in Asia (Rödenbeck et al., 2003a;Peters et al., 2007). Published inversion results also frequently disagree with estimates of carbon fluxes from biogeochemical models (Janssens et al., 2003;Peters et al., 2007;Potter et al., 2007). We will examine our results in this light, but concentrate on the uncertainty improvement of the added measurement sites in this paper.
In Section 2, we describe the estimation method. In Section 3, we present global and North American results for two typical global measurement networks and a third network including five additional continental sites. The results are discussed in Section 4. We conclude in Section 5 with recommendations for applicability of the method to future experiments.

Estimation method
The Bayesian synthesis inversion method used is shown in Fig. 1 and described by Enting (2002) and Tarantola (2005), with the experimental protocol (Gurney et al., 2000) following closely that of the TransCom interannual variability model intercomparisons (Baker et al., 2006;Gurney et al., 2008). We depart from the TransCom inversion method in a few important respects. (1) For observation data, we use monthly means and standard deviations derived directly from site observations of carbon dioxide (Rödenbeck et al., 2003b) instead of using a smoothed data product, such as GLOBALVIEW-CO2 (GLOBALVIEW-CO2, 2007) and assigned uncertainties (Baker et al., 2006). (2) We use annually varying meteorological driver data in our transport modelling.
(3) Finally, we include biomass burning emissions explicitly. The solution is for monthly carbon source/sink estimates (2000)(2001)(2002)(2003)(2004) for 47 subcontinental regions and ocean basins us-ing monthly mean carbon dioxide mixing ratio measurements. Figure 2 shows the region definitions and locations of the observing sites. The problem is ill-constrained due to the sparse and uneven distribution of observations; the method solves for adjustments to natural land-and ocean-atmosphere exchanges (referred to here as background fluxes) within the constraint of imposed prior uncertainties. Fossil fuel emissions and biomass burning emissions are assumed to be correct and not adjusted in the inversion process.
Following Baker et al. (2006), the atmospheric carbon dioxide at a measuring site can be represented as the linear combination of responses at the location to the background fluxes and to the unknown adjustment fluxes from each of the regions and months, where c obs is the time series of monthly carbon dioxide observations, c fwd is the modelled concentration time series using the background fluxes, H is a transport matrix (described below) and x are the unknown monthly adjustments to the background terrestrial and ocean fluxes for each of the 47 regions. Solving for the unknown monthly adjustments x is done, using a singular value decomposition approach for numerical stability (Rayner  , 1999), to the minimization of the cost function where R is the covariance matrix specifying the observation, transport and representation error of the data, x 0 are a priori estimates of the solution and P 0 is the covariance matrix of uncertainties of these a priori estimates. In this experiment both R and P 0 are described by diagonal matrices, as has been common practice in previous TransCom experiments Gurney et al., 2004;Baker et al., 2006). This cost function minimization is a least squares solution weighted by the variability of the observations and penalized for deviation from the a priori estimates. The analytical solution for the flux adjustmentsx iŝ and the a posteriori covariance matrix is given by allowing examination of the covariances for independence of the solution. This analytical a posteriori uncertainty is a function of the uncertainty of the prior flux estimate and the uncertainty attributed to the observations and assumes that the problem meets the requirements of the method. We have oversimplified both the data errors (uncertainty of observations and transport) and the model errors (prior flux uncertainty) by assuming they are independent and can be modelled using Gaussian distributions. Tarantola (2005) warns that the a posteriori solution may not be robust and may be particularly sensitive to outlier data points if these assumptions are not met. As a consequence, the accuracy of the a posteriori covariance is limited by the assumptions inherent in the a priori flux uncertainty and data uncertainty assignments. Even if the a posteriori covariance is optimistic, analysis of the region-region covariances provides valuable information about the flux estimates and the extent of the constraint provided by the network. Another product of this inversion method is a set of predicted observations, the carbon dioxide concentration values that would be expected at each observation site given the a posteriori flux solution. The predicted observations provide a basis for analysis of the data residuals, which cannot be done with the fluxes. The analytical solution delivers an adjustment to the background fluxes for each region and month in 2000-2004. Results for the central years 2001-2003 are retained for analysis. The first and final year of the solution are discarded to minimize edge effects, including inaccuracies introduced by the spinup method described below and the length of time it takes an atmospheric signal to reach distant observation sites. All results reported here are the adjusted land and ocean fluxes including the assumed biomass burning emissions for each month in order to compare with published results. The assumed background fossil emissions are not included. Choices for the components of the method shown in Fig. 1 follow.

Background fluxes
The terrestrial background flux is the annually varying hourly SiB3 terrestrial flux  as prepared for the TransCom continuous experiment Patra et al., 2008). An alternative terrestrial flux with monthly (but not diurnal or annual) variability, the CASA climatology  used in previous TransCom experiments, is also tested to investigate the solution dependence on the terrestrial background flux used and the time resolution of the flux. The ocean background flux is the monthly climatology of Takahashi et al. (2002), also as used in the TransCom continuous experiment. The seasonally varying fossil fuel emissions, as described by Erickson et al. (2008), are based on the 1995 spatial distribution of Brenkert (1998) and Li (1996), with annual totals scaled to the appropriate years based on Marland et al. (2007). The interannually varying monthly biomass burning emissions are from the Global Fire Emissions Database version 2 (GFED2) . Each background flux (4 tracers × 5 yr) is regridded to the horizontal spatial grid of the transport model and used as a surface boundary condition for a forward run from the beginning of the applicable year to the end of 2004. The model output is sampled hourly to allow for matching of sampling of the transport model output with the observations (see Section 2.8) when constructing the c fwd time series represented in eq. (1).

Transport response functions
The 47 regions specified for the solution are shown in Fig. 2 and listed in Table 1. The 11 ocean regions are the same as those used in the TransCom experiments. The 36 land regions are conformable to the TransCom land regions. The Boreal Asia, Temperate Asia, Tropical Asia and Australia regions are defined as in TransCom. Remaining TransCom regions are subdivided into subregions: Boreal North America (3 regions); Temperate North America (7); Tropical America (3), Temperate South America (2), Europe (7), Northern Africa (5) and Southern Africa (5). These continental regions are subdivided, roughly according to biome, to test the degree to which the solutions can be constrained given current and future observing sites. For North America and Europe the number of regions is similar to the number of available continental observing sites. The results for these regions are aggregated to larger regions for reporting and comparison with published results.
To construct the transport matrix H, emissions of known magnitude (1 month duration at an annual rate of 1 Pg carbon emission) are run forward as inert species through the transport model for each region-month (47 regions × 12 months × 5 yr), be-ginning in the applicable month and year and ending after 25 months of transport or the end of 2004. Model output is sampled hourly at each observation location through the 25 months (or less) of each region-month transport run, and then assumed to remain at a constant well-mixed, residual level after the transport is terminated through the end of 2004. Longer transport runs showed that the assumption of a well-mixed concentration beyond 25 months is justified.
The flux spatial distributions (or flux patterns) for the transport response functions within each terrestrial region are based on annual NPP simulated by CASA monthly climatology , scaled to the transport model grid cells within each region. The pattern is derived from the similar regional patterns in the TransCom experiments (Gurney et al., 2000) and regridded to the spatial resolution of the transport model used here. There is no variation by model grid cell within the ocean regions. These land and ocean flux patterns constitute a hard constraint in this inversion method. The inversion solution specifies the adjustment to the background flux for the region as a whole; the distribution of the flux within the region is fixed by these flux patterns.

Tracer transport model
The tracer transport model used in this experiment is the Colorado State University version of the NASA Parametrized Chemistry Tracer Model (PCTM), described in Kawa et al. (2004) which has participated in TransCom experiments. Winds, temperatures, diffusion coefficients and convective mass fluxes are from the NASA Goddard Earth Observation System 4 (GEOS-4) data assimilation system (Bloom et al., 2005). The 6-hourly meteorological driver data are linearly interpolated to the 15-min time step of the model. The model is run on a 2.5 • longitude by 2.0 • latitude grid, with 25 hybrid vertical layers. The model exhibited intermediate performance in the TransCom Interannual Variability model intercomparisons (Baker et al., 2006;Gurney et al., 2008) and is also one of the three models used in the  study of the potential bias in inversions caused by using non-varying background fossil emissions.  chose PCTM for its intermediate performance among the TransCom models with regard to key transport characteristics: the simulated annual zonal mean surface carbon dioxide concentration responses to non-varying fossil fuel emissions and, separately, to the neutral biosphere flux. The responses to these surface fluxes are dependent on the volume in which the surface fluxes are mixed (the planetary boundary layer depth), the seasonal phasing of the surface fluxes and the transport and the subgrid-scale parametrization used for vertical mixing. For PCTM, the interhemispheric gradient established by the fossil emissions is relatively large compared to the other TransCom models , but there is a relatively weak seasonal rectifier (covariance of the terrestrial background flux and seasonal differences in vertical mixing close to the model surface, resulting in an elevated annual mean mixing ratio in the atmospheric boundary layer relative to the free troposphere (Denning et al., 1995)). The interhemispheric exchange time, as defined in Denning et al. (1999), is slightly faster than average among the TransCom models (Kawa et al., 2004). Additional detail about the transport model performance can be found in Parazoo et al. (2008).
To establish the atmospheric interhemispheric gradient of carbon dioxide, we used a spin up procedure that differs from the technique used in TransCom-IAV (Baker et al., 2006), where there is no annual variability in the modelled responses used to construct the transport matrix H. The background fluxes were emitted as boundary conditions for 1 yr and the transport was continued for an additional 24 months for each model that participated in TransCom-IAV. In a similar manner the tracer response functions are the result of forward runs with one month of emission and the following 36 months of transport only. Each participating modeller made his own choice of transport fields. In constructing the transport matrix H, the same 36 months of responses are reused, beginning in each year of the inversion period. In TransCom-IAV, with its long analysis period, the first years are discarded in the analysis as they contain only the first years of responses. In our case, our goal was to use transport fields to match the years of the carbon dioxide observations. We had available only 5 yr of transport fields (2000)(2001)(2002)(2003)(2004) at the time the forward runs were executed. In order to maximize an already short analysis period, we augmented the year 2000 responses with representative model response data to simulate responses to background fluxes and monthly response functions from years before 2000. We then discarded inversion results from the year 2000 in our analysis results to minimize any potential effects of this procedure.

A priori constraints
The a priori constraint consists of a flux term and an uncertainty specification at monthly resolution for each region-month flux adjustment in the inversion solution and for the background fluxes. The region-month a priori flux adjustments in this experiment are set at zero, implying 'no correction' to the monthly background terrestrial and ocean fluxes, as well as the fossil fuel and biomass burning emissions. This procedure differs from some other methods that may incorporate best guess adjustments or biomass burning emissions or other land use change fluxes in the a priori specifications. The terrestrial background flux is annually neutral; the ocean background flux specifies a global annual sink of 1.6 Pg C; biomass burning emissions are ∼2 Pg C annually; fossil flux emissions are ∼7 Pg C annually for this period. In order to balance the global growth in atmospheric carbon dioxide concentration, there is an expectation that the inversion-adjusted terrestrial and ocean region-month fluxes will be a net sink. The region-month uncertainties used to populate the variances in the diagonal values of P 0 , the a priori flux covariance matrix, vary by month, summing globally to 5.4 Pg C annually. This range is larger than the 2.8 Pg C annual uncertainty used in the TransCom interannual variability control (Baker et al., 2006) and network sensitivity (Gurney et al., 2008) experiments. The magnitude of the monthly uncertainties is intended to be loose enough to allow the inversion to make substantial adjustments to the background terrestrial and ocean fluxes while remaining biogeochemically realistic. These region-month uncertainties are calculated as the sum of the magnitudes of three months of the applicable background flux (terrestrial or ocean) for the region centred on each month. This allows the inversion 'room' to correct for timing differences in seasonal cycles and generally allows for more latitude for adjustment in months when the fluxes are of largest magnitude.

Observation site networks
Observation site locations are shown on the map in Fig. 2. Table 2 details the observation site names, locations, responsible agencies, references and the range of monthly variability of observations. The observation sites include observatories, tall towers, flask sampling sites and an aircraft vertical Tellus 62B (2010), 5 and carbon dioxide time series from five flux towers. We have selected quasi-continuous measurement programs (in-situ measurements available as hourly averages) preferentially over discrete measurement programs (sampled approximately weekly in flasks for later analysis). We have excluded colocated measurement programs, retaining the highest samples from tall towers and the quasi-continuous measurement programs that are colocated with discrete measurement programs. The goal is to use measurements calibrated to the WMO standards for carbon dioxide (Zhao et al., 1997;Tans et al., 2003;Zhao and Tans, 2006). Calibration to WMO standards across multiple measurement methods and agencies has proved to be a challenge . Law et al. (2003) show that the impact of interagency calibration offsets can be accommodated in inversions on synoptic time scales (∼5 d), but may be significant on the monthly time scale we use in this experiment. To the extent that we have included observations from multiple measurement programs, we may have introduced bias into our results.
Stations are selected based on data availability to minimize bias introduced by gap filling. An upper limit of 12 missing months during the 2000-2004 time period was allowed, with a few exceptions for sites which began operation after 2000 (tall tower in Moody, Texas; flux towers at Southern Great Plains and Tapajos). This requirement prohibited the use of a number of sites that would be used if a smoothed, extrapolated data product such as GLOBALVIEW (GLOBALVIEW-CO2, 2007) were used. Three networks are tested in this experiment. The first is a 43-site base network, sourced from three agencies: National Oceanic and Atmospheric Administration Earth System Research Laboratory (NOAA ESRL), Environment Canada (EC) and CSIRO Marine and Atmospheric Research GASLAB in Australia. All of the observation sites in this network are also used in CarbonTracker (Peters et al., 2007;CarbonTracker 2008, http: //carbontracker.noaa.gov); we cannot use some of the CarbonTracker network sites due to the completeness of record requirement of our method. A 73-site enhanced network includes sites from other agencies, including high altitude aircraft and mountain top sites, archived at the World Data Center for Greenhouse Gases. These sites are used to improve global coverage. The last network, the 78-site continental extension network, adds the five flux tower observation sites. A primary objective of this experiment is to show the effect on the inversion estimates of including the four added North American continental sites.

Observation data preparation
Observations are selected with minimal screening for data quality issues identified by the observing agencies. For quasi-continuous continental sites, mid-day hours are selected (12-16 local standard time), except for mountain top sites where midnight hours are used. This selection is intended to maximize the contribution to the monthly mean of the hours representing wellmixed atmospheric conditions; these are conditions most likely to be modelled correctly in the tracer transport model. Monthly means are calculated from all available hours in the daily selection time periods without regard to meteorological conditions. Data uncertainties for each month are computed as the standard deviations of the available data for the month, with an imposed minimum of 0.5 ppm. These uncertainties, as variances, are used to populate the diagonal data covariance matrix R (eqs 2, 3 and 4). We have not included additional model-data mismatch error in R to account for errors in model transport or the representation of regions by point observation sites; we have equated synoptic monthly variability of the selected hours of observations with the total data error. Preliminary experiments with added model-data mismatch error effectively de-weighted observations to the point of overreliance on the a priori fluxes, with minimal adjustments to the background fluxes and large data residuals. This aspect of our method warrants further research. Table 2 includes the range of monthly variability for each station. Compared to the TransCom interannual variability experiment (Baker et al., 2006), these uncertainties are larger by ∼0.2 ppm at the minimum and up to two times larger at the maximum for remote sites. For continental and coastal sites, the minimum uncertainties are comparable and maximum uncertainties are at least as large as those in Baker et al. (2006), which used the GLOBALVIEW smoothed data product for the observation data time series and added model-data mismatch error for some observation sites.
For this analytical solution method, a regular time resolution is implied, making it necessary to fill gaps in the data for months with no observations. For observation sites that are represented in GLOBALVIEW (GLOBALVIEW-CO2, 2007), gaps are filled with the average of the GLOBALVIEW site data for the weeks corresponding to the month to be filled. For stations with no GLOBALVIEW representation, a climatology is constructed of the monthly departure of existing observations at the station from the GLOBALVIEW marine boundary layer (MBL) value for the station latitude. This difference climatology is applied to the MBL value for the missing months to fill the gap. The uncertainty for each gap-filled month is assigned as either the climatologic observation variability for the month for the station in 2000-2004 or a mean variability for all months with existing observations for the 2000-2004 time series for the station.

Model sampling
Model sampling is accomplished by saving vertical profiles of tracer concentration, pressure and temperature at the model grid cell nearest to each observing site hourly for the duration of Tellus 62B (2010), 5 the forward runs of the tracer transport model. Coastal stations are sampled in an off-shore grid cell to approximate observation protocol for flask sampling of background air. For surface observing stations the concentration samples are taken from the surface layer in the chosen grid cell. For elevated stations, noted as non-surface in Table 2, vertical pressure interpolation is used to resolve the concentration sample to an equivalent elevation. Hourly samples from all transport model forward runs were saved as described in Sections 2.2 and 2.3. The hourly samples that match the measurement time stamps of the observations are selected for use in the construction of the transport matrix H (eq. 2), a concept referred to as cosampling (Peters et al., 2007). For the discrete flask and aircraft observation sites, 3 h of model samples centred on the hour of each observation are used to smooth discrepancies in timing of weather events in the transport. Monthly mean model samples are calculated from these cosampled hours for each observation site from the forward runs for each of the region-month pulses and the background fluxes and used to populate the transport matrix H.
In the event of a month-long gap in the observation data, a default model sampling strategy is required. In this case model samples are chosen from the same hours that an observation sample would have been taken (daily selected hours for quasicontinuous sites and on 5 d during the month for flask sampling sites). The default hours for flask sampling are based on the distribution of hours when samples were taken during 2000-2004 at the respective sampling sites. Monthly means computed from these selected hours are used to fill in the gaps in the cosampled transport matrix H.

Results
All results shown are the adjusted land and ocean a posteriori flux estimates including the fixed biomass burning emissions, but not the fossil fuel emissions. This convention allows comparisons to previous work in the literature. Error bars in the figures are ± the mean annual uncertainty from the a posteriori covariance. CarbonTracker2008 results for the same time period are shown for reference (circles). All flux estimates include biomass burning emissions, but not fossil emissions. for this period. The proportion of the global sink attributed to land regions across the three networks (58% for the base network, 46% for the enhanced network and 53% for the continental extension) shows, at the global level, some sensitivity of the solution to the composition of the network. Although the tropical and southern regions are not individually well-constrained, for all three networks the tropical land/ocean source and the southern land/ocean sink sum to a 0.3 Pg C yr −1 source. The tropical plus southern latitudinal band is constrained by these networks, but the flux distribution within these bands is not well-defined. The ocean region fluxes are mainly consistent within uncertainties across the solutions obtained from the three networks. The three networks attribute 71-77% of the northern sink of 3.0-3.1 Pg C yr −1 to the land regions, with the sink for the ocean regions (North Pacific, North Atlantic and Northern Oceans) consistent at 0.7-0.9 Pg C yr −1 across the three solutions. However, the distribution of the land fluxes among continents depends very much on the network. For the base network, the land sink is split between Asia (Boreal and Temperate Asia regions accounting for 52%) and North America (Boreal and Temperate North America regions accounting for 37%), leaving Europe with 11%. Introducing east Asian observation sites, the western Pacific high altitude aircraft observations and mountain-top observations in Europe in the enhanced network shifts some of the Asian flux to Europe for an Asia:North America:Europe balance of 30%:39%:31%. This flux redistribution is accompanied by reductions in uncertainty for Europe and Asian regions (for example, Tropical Asia region annual uncertainty changes from 0.42 to 0.31 Pg C yr −1 and Europe from 0.54 to 0.48 Pg C yr −1 ). The continental extension network with the four additional continental North American sites shifts more of the flux to North America (22%:53%:25%), with a reduction in uncertainty in North America (from 0.51 to 0.40 Pg C yr −1 ), but not in Asia or Europe. The absence of observation sites in Boreal Asia and North Central Europe makes these regions under-constrained.

Global results
Two network-related shifts are seen in ocean fluxes between the base network and the enhanced network: the first is a transfer of source from the Tropical East Pacific and Tropical America to Temperate South America; the second is a shift of the distribution of the sink between the Southern Ocean and the South Indian Ocean with only modest reductions in uncertainty.   Southern Hemisphere observation site additions in the enhanced network include Jubany Bay in Antarctica, Cape Point in southern Africa and the southern branch of the high altitude western Pacific aircraft flights. The flux tower in the Amazon in the continental extension network makes little difference in either flux or uncertainty in the aggregated South American regions. The seasonality of the modeled atmospheric concentration of carbon dioxide at this site fits the data poorly before the inversion and the post-inversion predicted concentration here is the worst fit of all of the observation sites used. This issue of correctly modeling the terrestrial carbon flux in the Amazon is partially addressed by Baker et al. (2008). The poor model fit to the site observation and the lack of any other observations on the South American continent likely account for the lack of improvement in the certainty of estimated fluxes here.
A fully independent solution using this analytic method would have no off-diagonal entries in the a posteriori covariance matrix. Here we examine an annually summarized example to illustrate how well the continental network solution conforms to this expectation. The a posteriori covariance matrix for 2002 annual fluxes for the aggregated TransCom regions and the continental extension network is presented in Fig. 4. Error variances are reported on the diagonal. The off-diagonal values indicate the covariances between the annual flux solutions in different regions and, therefore, show the independence of the estimates. Shaded off-diagonal values, in variance units of (Pg C yr −1 ) 2 , impact the annual uncertainty of the aggregated regions by 0.1 Pg C or more; this cut-off is equivalent to the smallest annual error variance of the 22 aggregated regions. The largest covariance between ocean regions is slightly larger than -0.01 in units of variance between the Southern Ocean and the South Indian Ocean. The largest covariances are between Tropical America and Temperate South America (-0.37 in units of variance), North Africa and South Africa (-0.21), Temperate South America and South Africa (-0.14) and Europe and Boreal Asia (-0.12). The negative covariances among the Tropical America, Temperate South America, North Africa and South Africa aggregated regions, bounded by heavy black lines in Fig. 4, indicate that the estimates for the two continents should not be treated as independent. The total flux may be constrained, but the partitioning is not certain. Introducing the observation site in South America, where there were none in the enhanced network, is not enough to improve the overall continental constraint. In North America, however, the four added observations reduce the 2002 annual variance for Boreal North America from 0.16 to 0.09 (Pg C yr −1 ) 2 , for Temperate North America from 0.33 to 0.17   (Pg C yr −1 ) 2 and the covariance between them from −0.12 to −0.05 (Pg C yr −1 ) 2 (bounded by heavy black lines in the lower left corner of Fig. 4). Rödenbeck et al. (2003a) examined a posteriori flux errors from the perspective of the correlation coefficients of the long-term fluxes for a period in the late 1990s in their fig. 13, finding several of the same "pair-sum" relationships that we see here.  Fig. 6 for the base network (triangles), enhanced network (diamonds) and continental extension network (squares). The annualized prior uncertainty and mean annual posterior uncertainty for each region are also listed in Table 1. Figure 6 again demonstrates that uncertainty reduction is achieved in the regions local to the added observations, a finding consistent with results from grid-scale inversions (e.g. Peylin et al., 2005;Gourdji et al., 2008;Lauvaux et al., 2008;Mueller et al., 2008). In regions where there are no observation sites in any of these three networks (for example, the Subtropical region within Temperate North America and most regions in South America), there is no reduction in uncertainty. Introducing one or two sites in a previously unrepresented region can result in a significant reduction in uncertainty. The Northeast region within Temperate North America shows little improvement in the enhanced network compared to the base network (36% vs. 34%); adding two flux tower observation sites in the continental extension network results in an uncertainty reduction of 70% from the prior uncertainty with a mean annual posterior uncertainty of 0.25 Pg C yr −1 compared to 0.54 Pg C yr −1 for the enhanced network. The South American flux tower is on the border between the Northern and Southern Amazon regions within Tropical America. The percentage reduction in uncertainty in the Northern Amazon compared to the prior doubles in the continental extension network compared to the base and enhanced networks (16% for the base network, 18% for the enhanced  network, 36% for the continental extension). The Southern Amazon region shows no improvement at all. Parazoo et al. (2008), using the same transport model and analyzed meteorological fields, found that atmospheric mixing in this region is dominated by vertical convection rather than the horizontal transport by synoptic weather systems prevalent in the midlatitude continents. It is not surprising, therefore, that one observation site in South America constrains only the immediate region.

North American results
Regions within North America already partially constrained by local observations show some improvement, but less than for the first observation site. The Southwest, Central Plains and Eastern Boreal regions contain one observation site in the base network (Niwot Ridge, WKT tall tower and Fraserdale tower, respectively). Adding the aircraft vertical profile observations at Carr, Colorado, to the Southwest region in the enhanced network results in a modest improvement from 40% to 48% from the prior specification. After adding a flux tower to the Central Plains region and the Eastern Boreal region in the continental extension network, the percentage reduction of the uncertainty improves from 43 to 60% for the Central Plains and from 58 to 64% for the Eastern Boreal region.
Similar improvements in the covariances can be shown in the non-aggregated, 47-region 2002 annual covariance matrix for the continental extension network (not shown). The covariances among the regions within North America are reduced; the covariance between Eastern Boreal and Northeast Temperate regions is reduced from -0.06 to -0.02 (Pg C yr −1 ) 2 , and the variances for Central Plains and Northeast are improved [0.11-0.05 (Pg C yr −1 ) 2 for the Central Plains and 0.30-0.05 (Pg C yr −1 ) 2 for the northeast]. Covariances between the North American regions and non-North American regions are reduced except for a negative covariance between Western Boreal North America and Boreal Asia. Results are more modest in South America, where the only appreciable improvement is the reduction in the covariance between Northern and Southern Amazon regions from −0.11 to −0.07 (Pg C yr −1 ) 2 . Covariances within South America and between South American and African regions are reduced but not eliminated. We cannot justify the inversion solution for the subregions within South America as independent with the networks tested.
The monthly flux solutions for two Temperate North American regions for the three networks are shown in Fig. 7. The base network results are shown in blue, the enhanced network in cyan and the continental extension in red. Dotted lines indicate the extent of the constraint provided by the prior uncertainty. Shading represents the a posteriori uncertainty (from the diagonal of the posterior covariance matrix) of the solution for the continental extension network. On the left in Fig. 7 is the Pacific Northwest region which has no local observation sites in any of these networks, and is upstream of all the added North American   observation sites. Other than a correction to the timing of the seasonal cycle of the terrestrial background flux, there is little change from the background flux and no difference in the a posterior flux estimate for the three networks tested. The Northeast region is shown in the right panel of Fig. 7. This region has no local observations in the base or enhanced networks, and the two flux towers with well-calibrated CO 2 measurements in the continental extension network. The a posteriori flux estimate for the continental extension network is smoothed, of smaller amplitude and reduced uncertainty relative to the prior. Similar uncertainty reductions that can be attributed to the presence of local observations is evident in the monthly flux solutions for Central Plains and Eastern Boreal regions, where the two other North American flux towers in the continental extension network are located. Figure 8 shows  (Waple et al., 2003), but a strong sink anomaly in the wet summer of 2003 (Levinson and Waple, 2004). The variability in the Northeast region exhibits the reverse behavior for these 2 yr. The 47-region covariance matrix shows an annual covariance of −0.01 (Pg C yr −1 ) 2 between these two regions. The negative covariance, although small, might be an indication of some dipole behavior between these two regions. All of the anomalies are within 1 standard deviation of the 3-yr means.

Impact of added measurement sites
The added continental measurement sites contribute to local uncertainty reduction without undue disruption to the global flux solution. We choose to use mid-day observations and matching transport model samples for continental surface sites based on comparisons of mid-day and all-hours sampling for both the observations and transport model output. As these mid-day hours are most likely to represent the desired well-mixed background state, they are also likely to represent an area much broader than the immediate vicinity of the observation site. This subsampling avoids the issue of nocturnal boundary layer representation in the transport models and driver meteorological data. Representation error is still possible; choice of site is critical. The competing roles of multiple sites in the same solution region (e.g. Howland Forest and Harvard Forest in the Northeast region of Temperate North America) warrant further research. There are a number of additional North American flux towers where carbon dioxide measurements are now well-calibrated. These present an opportunity for improving the uncertainty of the North American carbon balance as determined by atmospheric inversions. The use of these stations is, however, dependent on the continuity of the measurement time series, and the degree to which the site is representative of the area for flux estimation. The inversion method used here requires continuous observation time series and an unchanging network. Other methods may be more amenable to the addition and deletion of observation sites for an inversion solution for long periods of time.

Results using an alternative background terrestrial flux
We used an alternative background terrestrial flux, the CASA monthly mean flux, to compare the a posteriori results for the continental extension network with the results using the SiB3 terrestrial flux which has diurnal, monthly and annual variability. The results for the 22 aggregated TransCom regions are documented in Table 3 and as squares in Fig. 9. Ocean region results are similar, as might be expected. The two background fluxes have differences in both amplitude and timing of the seasonal cycle in North America. In well-constrained regions, such as the Northeast Temperate North America, both background terrestrial fluxes are adjusted to a result that is consistent within the a posteriori uncertainty. The results for CASA and SiB3 are less consistent for under-constrained regions; Tropical America is a good example of this. The sensitivity to the background fluxes warrants further research.
We have found the model samples from the forward runs of SiB3 overestimate observed seasonal amplitudes of carbon dioxide concentrations by as much as 40% at observation sites in some northern temperate and boreal regions in North America. This can be seen in the prior uncertainty ranges (derived from the SiB3 fluxes) in Fig. 7, especially in the Northeast Temperate North America region. The CASA terrestrial background flux has smaller annual amplitude in these same regions, and the modelled seasonal amplitudes are within 10% of the observed amplitudes at the same locations that are overestimated by SiB3. As can be seen in Fig. 9, the inversion solutions using the alternative background fluxes are similar in regions that are well constrained.

Comparison to published results
The results of this inversion are placed in the context of other contemporary inversions for the same time period in Fig. 9 and Table 3. Included in the table and figure are the SiB3 and CASA versions of our experiment (squares in Fig. 9), Car-bonTracker 2008 (circles) from http://carbontracker.noaa.gov, model mean results (diamonds) from the 95-site network in the TransCom Interannual Variability (IAV) network sensitivity study (Gurney et al., 2008), and the NASA GSFC PCTM model results (triangles) from the same IAV study. All results are for [2001][2002][2003]. The uncertainties shown in Table 3 and Fig. 9 for the TransCom IAV model mean are the total error including the model spread. The 95-site network used in the TransCom IAV network sensitivity study is similar to our enhanced network, except that the GLOBALVIEW smoothed data product is used for observation data and colocated measurement programs are not excluded. In Fig. 9 we show inversion solutions using different data sources, different transport models and meteorological driver data, different background ocean and terrestrial fluxes, and different inversion methods. One very obvious difference in the partitioning results in Fig. 9 is in the greater magnitude of the uncertainties estimated by CarbonTracker. This is a methodological difference and an attempt to incorporate transport model error more comprehensively (Peters et al., 2007) when using only one transport model. Other studies have suggested that the transport model used and the uneven and sparse spatial distribution of observation sites are equal contributors to the differences in results Rödenbeck et al., 2003a,b;Gurney et al., 2008). The other a posteriori uncertainties shown in Fig. 9 are derived analytically (with the addition of model spread for the TransCom-IAV model mean). In general our SiB3 and CASA results, which share the same transport, network composition and data selection protocol, are much closer to each other than to the other inversions in this comparison. Figure 9 shows that the TransCom IAV-PCTM result tracks very closely to the TransCom IAV model mean except in the partitioning between North America and Extratropical Eurasia where the TransCom IAV-PCTM more closely resembles the results of our study. Our results in the partitioning between Global Land and Global Ocean are more like CarbonTracker than the TransCom IAV model mean and TransCom IAV-PCTM results with greater land sink than ocean sink. All of the inversion results report a net source for the combined Tropical and Southern Land.
In contrast to CarbonTracker, however, our results place more of the northern land sink in North America than in extratropical Eurasia (Europe, Boreal Asia and Temperate Asia) in agreement with the TransCom IAV model mean and TransCom IAV-PCTM. Only the CarbonTracker results in Fig

Adjustments of the background ocean fluxes
Our inversion results are generally consistent in zonal distribution with the ocean air-sea flux estimates of Takahashi et al. (2002), which we use as our background ocean flux. Our flux solution reduces the net global ocean sink from 1.6 Pg C to 1.2 Pg C annually, with smaller sinks in temperate waters and smaller sources in the tropics. Our experiment reduces the Southern Ocean sink specified by Takahashi et al. (2002), but not by as much as CarbonTracker (see Table 3) or the update to these airsea fluxes (Takahashi et al., 2009) for the climatologic year 2000, which shows little or no net annual flux in the Southern Ocean due to offsetting effects of biological drawdown and upwelling of carbon-rich waters.

Subcontinental flux partitioning
The relatively large (with respect to CarbonTracker) Temperate North America terrestrial sink found in this study is not evenly distributed across the continent. Table 1 and Fig. 5b) suggest that the northeastern U.S. is responsible for more than half of this sink. The additional data constraint provided by the flux tower mixing ratio time series suggests that, within the limits of the uncertainties assessed here, this finding is significant. It is possible that this could be the result of transport errors, network biases, or other systematic errors that this study is not able to resolve. It is also notable that in Europe, the other relatively densely instrumented continent, the North Central region presents a dominant sink and the Eastern region a source, both significant with respect to the bounds of our uncertainty estimates. These subcontinental results show the promise of enhanced continental networks and, while preliminary, deserve further study.

Uses of the posterior covariance matrix
We have used the analytic product of the inversion method, the a posteriori covariance matrix, in three different ways in this experiment. We have used the magnitudes of the matrix elements as an indication of the uncertainty of the flux estimates. We have also used a change in magnitude as an indication of the value of adding sites to the network. Finally, we have used the a posteriori covariance to assess the spatial correlations between regions in the flux solutions. As we noted in Section 2.1, the accuracy of the a posteriori covariance is limited by the assumptions inherent in the a priori flux uncertainty and model-data errors, both of which are difficult to test. The magnitudes of the a posteriori uncertainties, as reported in Tables 1 and 3 and in the figures, are dependent upon the accuracy of our assumed prior uncertainties. The other uses of the a posteriori covariance should, however, be largely independent of the magnitude of the prior uncertainty estimates and thus more robust. CarbonTracker (Peters et al., 2007) includes a broader uncertainty assessment including some analysis of uncertainty caused by atmospheric transport. This may explain the larger error bounds in the CarbonTracker flux estimates.

Conclusions
We have shown the value of additional continental observations in a standard Bayesian global atmospheric inversion. We have removed some of the methodological simplifications used in previous Bayesian synthesis inversions by using transport fields appropriate to the years of the observation data, by using real observation data, and by explicit inclusion of biomass burning emissions. We have also incorporated careful examination of the a posteriori covariance in the evaluation of the quality of the inversion output. Specific findings include the following: 1. Using a judicious combination of region size and location of additional continental observation sites, the a posteriori flux uncertainty and region-to-region covariance are significantly reduced.
2. Flux towers are typically sites in representative ecosystems and have existing measurement infrastructure. With some additional instrumentation and a careful calibration of carbon dioxide concentrations to global standards, these sites can be used in the global measurement network used for atmospheric inversions. With the extent of the current global flux network, there is a major opportunity for improved coverage in the observation networks used for global and regional atmospheric inversions.
3. The transport model and analysed meteorological fields appear to have more impact on our flux estimates than the added observation sites or the background terrestrial flux. For the well-constrained aggregated region of North America, the TransCom-IAV PCTM results are more like our findings than CarbonTracker. Transport uncertainty is an important problem that this study does not address. 4. We confirm the findings of other researchers (e.g. Baker et al., 2006;Jacobson et al., 2007a,b) that the tropical and southern land regions are underconstrained to the extent that the inversion source/sink estimates for these regions cannot be considered independently. The observation networks tested in this experiment do not improve this situation. An increase in observation density, for example from satellite observations, is required to change this situation. 5. In our experiment, we find that the Northeast subregion of Temperate North America dominates the North American sink. While this is a preliminary finding, based on our limited experiment, it does suggest that increasing both the coverage of the observation network and the number of regions in the inversion, can help to localize critical source and sink regions for the global carbon budget.