The covariation of Northern Hemisphere summertime CO 2 with surface temperature in boreal regions

. We observe signiﬁcant interannual variability in the strength of the seasonal cycle drawdown in northern midlatitudes from measurements of CO 2 made by the Total Carbon Column Observing Network (TCCON) and the Greenhouse Gases Observing Satellite (GOSAT). This variability correlates with surface temperature in the boreal regions. Using TCCON measurements, we ﬁnd that the slope of the relationship between the X CO 2 seasonal cycle minima and boreal surface temperature is 1.2 ± 0.7 ppm K − 1 . Assim-ilations from CarbonTracker 2011 and CO 2 simulations using the Simple Biosphere exchange Model (SiB) transported by GEOS-Chem underestimate


Introduction
Fossil fuel burning, the oceans, and the global terrestrial biosphere control the atmospheric concentrations and variability of carbon dioxide (CO 2 ). On interannual timescales, variations in fluxes from terrestrial ecosystems, driven by changes in surface temperature and precipitation, are the dominant drivers of variability in atmospheric CO 2 (Francey et al., 1995;Langenfelds et al., 2002;Rayner et al., 1999;Welp et al., 2011;Cox et al., 2013). On these timescales, warm years tend to be associated with more rapid increases in atmospheric CO 2 , and cool years with reduced growth rates (Braswell et al., 1997). This positive relationship between temperature and atmospheric CO 2 is attributed primarily to an increase in ecosystem respiration (R e ) with increasing surface temperature, and a concurrently muted gross primary production (GPP, or photosynthesis) (Doughty and Goulden, 2008).
The shape of the seasonal cycle in atmospheric CO 2 in the northern midlatitudes is primarily determined by the seasonal imbalance between R e and GPP (Randerson et al., 1997;Messerschmidt et al., 2013), which is referred to as net ecosystem exchange (NEE), where positive NEE is defined here as a net flux of CO 2 into the atmosphere. At high latitudes, there is a significant time lag between GPP and R e : GPP peaks around the summer solstice when photosynthesis is at a seasonal maximum, whereas R e peaks later in summer when air and ground temperatures are warmest. This creates a negative NEE during the growing season (June, July and August) (Lloyd and Taylor, 1994). The growing season NEE has the largest magnitude over the temperate and boreal forest region of the Northern Hemisphere , producing the observed seasonal cycle in northern midlatitude CO 2 (Machta, 1972, and Fig. 1). The seasonal cycles at all northern midlatitudes are highly sensitive to changes in the NEE in the boreal forests (D'Arrigo et al., 1987;Randerson et al., 1997;Keppel-Aleks et al., 2011), where there is also significant temperature-driven variability (Randerson et al., 1999;Piao et al., 2008). Thus, the atmospheric seasonal cycle and its variability can provide insight into the response of carbon pools to climate (Keeling et al., 1996;Randerson et al., 1997).
Previous analyses of the role of temperature on atmospheric CO 2 have relied on highly precise and accurate atmospheric CO 2 concentrations or fluxes measured by surface or near-surface in situ instruments located throughout the world. Tower flux measurements have shown a strong sensitivity of carbon exchange to surface temperature at several locations, including Niwot Ridge, Colorado (Sacks et al., 2007), Howland Forest, Maine (Richardson et al., 2007), and in the North American prairie region (Arnone et al., 2008). In this paper, we use measurements of column-averaged dry-air mole fractions of CO 2 , denoted X CO 2 , from the Total Carbon Column Observing Network (TCCON, Wunch et al., 2011a) and from the Greenhouse Gases Observing Satellite (GOSAT, Hamazaki, 2005;Yokota et al., 2009) to examine the interannual variability of the seasonal cycle minimum and its relationship with surface temperature. Compared with surface in situ measurements of CO 2 mole fractions, X CO 2 is influenced much less by planetary boundary layer height changes, and possesses a much larger spatial sensitivity footprint (on the order of hundreds to thousands of kilometers, Keppel-Aleks et al., 2011). The north-south gradients in X CO 2 in the Northern Hemisphere summertime are large, and hence the latitudinal origin of the measured air parcel strongly influences its X CO 2 (Keppel-Aleks et al., 2012).
There have been marked interannual differences in the seasonal cycle minima of X CO 2 in the Northern Hemisphere in recent years (Figs. 1, 2). Guerlet et al. (2013) have also described significant differences in the 2009 and 2010 seasonal cycle amplitudes using the GOSAT data. These interannual differences are sometimes underestimated in Earth system models, which explicitly represent feedbacks between climate and terrestrial carbon fluxes (Keppel-Aleks et al., 2013).
We find that the measured X CO 2 seasonal cycle minima are positively correlated with the measured surface temperature anomalies in boreal regions. Using atmospheric total column measurements, the CarbonTracker assimilation output (CarbonTracker, 2011) and global CO 2 simulations using the Simple Biosphere model (Sellers et al., 1996a), we investigate the following processes for their contribution to the observed interannual variability in X CO 2 : the relative contributions of fire, fossil fuel, terrestrial biosphere, and ocean fluxes, and dynamical drivers of interannual variability.
In the following sections, we describe the data and model sources, and outline our analysis methods. We then discuss the results of the analyses.

TCCON
The TCCON is composed of ground-based Fourier transform spectrometers distributed throughout the world that provide measurements of X CO 2 (Wunch et al., 2011a). We use data from the four longest-running Northern Hemisphere TCCON sites: Park Falls, Wisconsin, USA (46 • N, 90 • W, Washenfelder et al., 2006); Lamont, Oklahoma, USA (37 • N, and Bremen, Germany (53 • N, 9 • E). In Bremen, the construction of the solar tracker was completed in 2006, so we include measurements from the beginning of the subsequent calendar year. The X CO 2 values from these four sites have been tied to the World Meteorological Office scale through comparisons with aircraft profiles Wunch et al., 2010;Messerschmidt et al., 2011). We use the GGG2012 version of the TCCON data, available from http://tccon.ipac.caltech.edu, applying the corrections to Park Falls (+0.8 ppm after 23 June 2011) and Lamont (−0.5 ppm after 14 April 2011) recommended on the TC-CON data description webpage: https://tccon-wiki.caltech. edu/Network_Policy/Data_Use_Policy/Data_Description.

GOSAT
The GOSAT satellite, carrying the Thermal and Near-Infrared Sensor for carbon Observation Fourier transform spectrometer (TANSO-FTS), was launched in January 2009 and has a ground-repeat cycle of 3 days and a footprint size of approximately 100 km 2 (Yokota et al., 2009;Crisp et al., 2012). We use the X CO 2 derived from GOSAT spectra by the Atmospheric CO 2 Observations from Space (ACOS) build 3.3 retrieval algorithm (Crisp et al., 2012;O'Dell et al., 2012). Data from 1 April 2009 through 13 April 2013 are used in this study. The data were filtered and corrected for retrieval biases using the method described in the ACOS Level 2 Standard Product Data Users Guide, v3.3, available from http://disc.sci.gsfc.nasa.gov/datareleases/acos-version-3.3.

D. Wunch et al.: Covariation of X CO 2 with boreal temperature 9449
We define GOSAT-TCCON coincidences in the same manner as Wunch et al. (2011b). We use relatively wide latitude (±5 • ), longitude (±30 • ) and time criteria (±5 days), but restrict the GOSAT measurements to those having a freetropospheric temperature within ±2 K of that measured over the TCCON station. This allows averaging of measurements of air with similar dynamical origin. All GOSAT data that satisfy these criteria for a given day are averaged.

CarbonTracker
CarbonTracker release 2011 (Peters et al., 2007, http:// carbontracker.noaa.gov/; henceforth CT2011) is an ensemble data assimilation scheme that uses surface, tower, and shipborne in situ measurements of atmospheric CO 2 and the TM5 atmospheric transport model to produce 4-D fields of CO 2 (CarbonTracker, 2011). The TM5 model (Krol et al., 2005) is driven by the European Centre for Medium-range Weather Forecast (ECMWF) assimilated winds (Uppala et al., 2005;Molteni et al., 1996). There are four CO 2 flux "modules" embedded within CT2011: one for each of fire, fossil fuel, terrestrial biosphere, and ocean. These fluxes add to a background field to produce variability in assimilated CO 2 . There are two priors for each of the fossil fuel, biosphere, and ocean flux modules, resulting in eight separate inversions, and the ensemble mean is reported as the result.
Only the terrestrial biosphere and ocean modules are optimized in the assimilation scheme: as with most assimilation schemes, the fire and fossil fuel fields are prescribed. The fire emissions are prescribed using the Global Fire Emissions Database (GFEDv3, CarbonTracker, 2011;Giglio et al., 2006;van der Werf et al., 2010;Mu et al., 2011). CT2011 assumes that the fossil fuel emissions are known from reported annual national and global inventories compiled by the Carbon Dioxide Information Analysis Center (CDIAC, Boden et al., 2011;CarbonTracker, 2011).
The terrestrial biosphere module in CT2011 is optimized based on priors from the monthly mean Carnegie-Ames-Stanford Approach (CASA) ecosystem exchange model (Potter et al., 1993;Randerson et al., 1997), which uses satellite measurements of the normalized difference vegetation index (NDVI) and fractional photosynthetically active radiation (fPAR) as proxies for plant phenology, and yearspecific weather. Diurnal and synoptic variability in R e is imposed through a Q 10 relationship with surface air temperatures (R e ∝ Q (T −T 0 )/10 10 ), assuming a Q 10 of 1.5 for respiration globally (CarbonTracker, 2011). The terrestrial biosphere fluxes and fire emissions described above were generated from the same versions of the CASA GFED model.
The ocean module optimizes prior fluxes provided by oceanic flux inversions , and measurements of partial pressure CO 2 in the ocean surface (Takahashi et al., 2009). The prior fluxes have a smooth trend, but no interannual variability. Any interannual variability in the optimized fluxes from the oceanic module of CT2011 is due to atmospheric surface winds interacting with the ocean surface, affecting the gas transfer efficiency. The magnitude of the interannual variability is small compared with the biosphere (CarbonTracker, 2011).
We use CT2011 output sampled at the locations of the TC-CON stations in this study. We smooth the CO 2 profiles with the TCCON column averaging kernels and a priori profiles, using the Rodgers and Connor (2003) method, and integrate the smoothed profiles to produce daily CTX CO 2 . (Note that to compute total columns from the CO 2 profiles associated with the individual flux modules, no smoothing is applied: we simply integrate the CO 2 profiles to produce CTX module CO 2 .) This version of CarbonTracker has recently been replaced because of the discovery of a bug in the TM5 transport model. This bug is thought to have a negligible affect on the CO 2 mole fractions and no effect on fluxes from fossil fuel and wildfire emissions (CarbonTracker, 2011).

GEOS-Chem and SiB
GEOS-Chem is a global chemical transport model with CO 2 simulations described by Suntharalingam et al. (2003Suntharalingam et al. ( , 2004 and updated by Nassar et al. (2010). Typically, the CASA ecosystem exchange model is used to generate a priori biospheric CO 2 fluxes in GEOS-Chem. Messerschmidt et al. (2013) have recently shown that replacing CASA with the Simple Biosphere model (SiB3, Baker et al., 2008;Sellers et al., 1996a;Parazoo et al., 2008) significantly improves the CO 2 seasonal cycle amplitude and phase compared with TCCON observations. SiB3 calculates year-dependent fluxes using satellite measurements of plant phenology (Sellers et al., 1996b) yielding significant interannual variability in GPP (Baker et al., 2010). Here, phenological parameters are prescribed from the Moderate Resolution Imaging Spectroradiometer (MODIS, Zhao et al., 2006). Ecosystem respiration in SiB is driven by a Q 10 relationship with surface air temperature (T ), modified by a soil moisture term, g(m) (Denning et al., 1996): The soil moisture term is defined in Denning et al. (1996, Eq. 8) and is related to the fraction of root zone soil porosity holding water, prescribed for each soil type by Raich et al. (1991). Based on the work of Raich and Schlesinger (1992), Q 10 in SiB is set to 2.4. Here, the GEOS-Chem model was run twice: once with SiB 2009 fluxes (referred to as "SiB 2009"), and once with year-dependent SiB fluxes (referred to simply as "SiB"). CO 2 profiles from the model are interpolated to the locations of the TCCON stations, smoothed with the TCCON column averaging kernels and a priori profiles, integrated, and then averaged to produce daily SiBX CO 2 .

Methods
To determine the seasonal cycle minimum date and value, we fit the measured X CO 2 and the CTX CO 2 time series with an annual periodic function superimposed on the Earth System Research Laboratory (ESRL) global annual CO 2 growth rate (Table 1, Conway and Tans, 2013). We do not assume the ESRL global annual growth rates for the GEOS-Chem time series; instead, we add an additional linear increase term (αx) to the fit. The fitted curves do not permit inter-annual variations in the seasonal cycle shape. The functional form of the fitted curve is a Fourier series:  Figure 2 shows the time series detrended using the ESRL global annual growth rates, marking the date of the seasonal cycle minimum with symbols.
To assess how well the models compare with the TCCON data, we follow Yang et al. (2007) and use a least squares fit to find the amplitude and phase that minimizes differences between the seasonal cycles computed from the models and the TCCON data. Figure 3 shows the detrended seasonal cycles from the fitted curves (i.e., the k = 1, 2 terms from Eq. 2) to the models and the data. The agreement is generally good: the amplitudes of the SiBX CO 2 seasonal cycles are too large  Table 3 for the amplitude ratios and time lags.
by ∼10-20 %, but show good agreement in the timing of the seasonal cycle (within 7 days). The amplitudes and time lags from the CTX CO 2 seasonal cycles match the TCCON data well, with the exception of an 8-day time lag at Lamont ( Table 3). The seasonal cycle minimum value (which we will call the "drawdown value") is determined by subtracting the fitted curves from the time series, and averaging the resulting X CO 2 within ±15 days of the seasonal cycle minimum. The reported errors on the drawdown value represent the standard deviation (1σ ) of the measurements. In years with fewer than 3 days of measurements near the seasonal cycle minimum, this averaging date range is extended to be within ±25 days (Białystok and Bremen data in 2010). The errors for these data points are tripled to reflect the additional uncertainty. (2). The parameters α and a 0 represent the linear secular increase: α is in ppm yr −1 , and a 0 is the y intercept in ppm. Parameter a 1 (ppm) represents the amplitude of the cos(2πx) term, and b 1 (ppm) represents the amplitude of the sin(2π x) term. The a 2 and b 2 parameters are the amplitudes of the cos(4πx) and sin(4πx) terms, respectively. The drawdown date (DD) and seasonal cycle maximum date (MD) are the dates of the local minima and maxima of the fitted curves, respectively, in day of the year. The large negative values of a 0 for SiB and SiB2009 are the result of using a straight-line fit to the data (i.e., α) instead of subtracting the ESRL growth curve, and fitting the detrended curve.
Site α a 0 a 1 a 2 b 1 b 2 DD MD Park Falls (TCCON) ESRL (see Table 1  To investigate the impact of the individual CT2011 modules (fires, fossil fuels, terrestrial biosphere, ocean) on the CTX CO 2 interannual variability, the fitting method described above (Eq. (2)) is not used, because the ESRL growth rate is only applicable to the total CO 2 , and there may be no periodicity to some of the individual components. The CTX module CO 2 anomalies are instead computed by subtracting a yearly CTX module CO 2 mean that has been linearly interpolated to each time step. The drawdown values are the averages of these anomalies within ±15 days of the seasonal cycle minimum  Figure 4 shows the 2009 and 2010 August temperature anomalies as examples. As we wish to evaluate the coupling between the X CO 2 and the CO 2 fluxes, we weight the year-specific August surface temperature anomaly by the 2009 integrated growing season respiration from SiB (R 2009GS e , in kg C m −2 ) between 30 • N and 60 • N. These values are then divided by the 2009 integrated growing season (June, July, and August) ecosystem respiration between 30 • N and 60 • N to compute a respirationweighted temperature anomaly, δT , for each year y. This weights the temperature anomalies more strongly in locations where the biosphere is active, and de-weights regions in which the biosphere is less active (i.e., over ocean, barren, or snow-covered areas).
where i is the latitude, j the longitude, and a ij the grid area (in m 2 ). The value δT has units of temperature (K). The correlation between the drawdown value and respiration-weighted temperature anomaly for Park Falls is plotted in Fig. 5, and the slopes ( ∂ X CO 2 ∂δT ) for all TCCON sites are listed in Table 4. The slopes for Park Falls, Białystok, Bremen and Lamont are consistent within their standard errors, giving a weighted average of 1.2 ± 0.7 ppm K −1 . The Bremen slope is significantly smaller than the other slopes, possibly because of Bremen's proximity to urban fossil fuel emissions. Excluding Bremen from the weighted mean results in an average of 1.4±0.8 ppm K −1 . The Białystok slope is consistent with the Park Falls and Lamont slopes, but due to low data yields during the summers of 2010 and 2012, the error in the estimated slope is large. The standard errors are generally large at all sites because there is significant day-to-day variability in X CO 2 near the seasonal cycle minimum due to the influence of synoptic-scale activity on the measured X CO 2 (Keppel-Aleks et al., 2012). Figure 6 illustrates this increased summertime synoptic-scale "noise".  Fig. 6. Detrended X CO 2 seasonal cycles measured at the four TC-CON stations. The colors are the free-troposphere temperatures (at 700 hPa, denoted T700). The color bar is scaled to emphasize the spring and summer months. The summertime X CO 2 is lowest when the T700 values are coolest, and highest when the T700 values are warmest. The inset is Park Falls data alone, illustrating the increase in synoptic-scale "noise" in the summer.  Fig. 7. A map of the mean X CO 2 −δT slope over the 30-70 • N region from two summertime drawdown periods using GOSAT data. The TCCON stations involved in this study are indicated with black stars. The panel on the left shows the zonal-mean X CO 2 −δT slope as a function of latitude. Note that the latitude scales are not the same for the left panel and map: the left panel is linear in latitude, while the map is a Miller projection.

Results and discussion
In the summertime, air that has originated from the north (cooler air) tends to have lower X CO 2 values, whereas air that originated from the south (warmer air) has higher X CO 2 values, giving rise to high variability in X CO 2 in the growing season. Covariations between the drawdown values and 30 • N to 60 • N surface temperatures (not weighted by R e ) show larger relative errors on the slopes, but consistent slope values (within error). Using the Mauna Loa growth rate instead of the global mean ESRL growth rate decreases the overall slopes by ∼ 0.2 ppm K −1 . The GOSAT slopes (0.7 ± 0.8 ppm K −1 ) are consistent with the TCCON values within the errors. More revealing from the GOSAT analysis is the spatial distribution of the GOSAT slopes and their latitudinal gradient, which decreases from ∼ 1.3 ppm K −1 near 30 • N to ∼ 0.7 ppm K −1 by 60 • N (Fig. 7). The spatial pattern is mostly uniform, except near the oceans, and in regions of large local fossil fuel emissions (e.g., western Europe). This is consistent with our understanding that X CO 2 has a very large spatial footprint that is essentially hemispheric on seasonal timescales. The cause of the spatial variability near the coasts is unclear.
The significant correlation of X CO 2 with temperature could point to a large-scale dynamical effect, fires, fossil fuel use, or a biospheric reaction to the temperature changes. It is unlikely to be related to oceanic fluxes, as the interannual variability in CO 2 from ocean fluxes is negligible ( Table 5). The possible contributing effects and an estimate of their relative importance are discussed in turn below.

Dynamics
Variability in the dynamical mixing of CO 2 within and between the midlatitudes and the tropics is expected to be correlated with surface temperature changes, since the meridional thermal structure both responds to and drives north-south transport of air (Trenberth, 1990;Webster, 1981). Therefore, it is plausible that the slopes seen in the TCCON and GOSAT data are influenced by interannual variability in the atmospheric mixing. To test this, we use the GEOS-Chem dynamical fields with static SiB 2009 biospheric fluxes, holding the NEE cycle constant from year to year. The results are shown in Fig. 5 for Park Falls, and in Table 4 for all sites. The SiB run with 2009 fluxes shows a weakly positive slope that is approximately 40 % of that observed, indicating that variability in the mixing does contribute to the observed interannual variability in X CO 2 minima. Because the covariations are computed with column-averaged CO 2 , which is less influenced by local dynamics, this suggests that large-scale dynamics are important: warm years are associated with enhanced mixing of high CO 2 air from the subtropics with low CO 2 air from the Arctic, and weaker mixing in cool years.

Fires
High temperature (and lower humidity) conditions might be expected to be correlated with wild fire activity (e.g., Westerling et al., 2006) and therefore increased atmospheric CO 2 (Zhao and Running, 2010;Patra et al., 2005). To determine whether variations in fires have a significant effect on the δT (K) Fig. 8. The top panel shows the TCCON X CO 2 calculated from measured X CO during the summertime drawdown against the respiration-weighted surface temperature, to probe the effects of forest fires on the X CO 2 . The solid lines show the best-fit lines, which have very poor fit statistics (R 2 < 0.2). The middle panel shows the CT2011 fire contribution to the summertime drawdown in CTX CO 2 , which also has a negligible relationship with the respiration-weighted surface temperature. The bottom panel shows the CT2011 fossil fuel contribution to the summertime drawdown, which has a negligible relationship with the respiration-weighted surface temperature.
variations in X CO 2 seasonal cycle minima, we analyze TC-CON measurements of X CO , a fire tracer measured simultaneously with the X CO 2 . The anomalies in X CO are calculated in an identical manner to X CO 2 , except that we do not use the ESRL CO 2 mean growth rate, but fit an additional linear increase parameter (αx). Following Akagi et al. (2011), we assume that the modified combustion efficiency of producing CO 2 from biomass burning is ∼ 88 % in the boreal forest, and we convert the measured X CO at the time of the X CO 2 drawdown into the X CO 2 contribution from fires. We show the relationship between X CO 2 caused by fires with the respiration-weighted temperature anomaly in Fig. 8. The variability in X CO 2 caused by fires is at most ∼ 0.05 ppm K −1 . Although CO is not a conserved tracer (its oxidation by OH leads to a lifetime of ∼ 1-2 months (Yurganov et al., 2004)), the small slope, even if a lower bound, suggests that fire does not contribute significantly to the observed ∼ 1.2 ppm K −1 variability. Consistently, the CT2011 fire signature (CTX fires CO 2 ) accounts for only 6-13 % of the total interannual variability in summertime drawdown CTX CO 2 in CarbonTracker (see Table 5). The fire anomalies ( CTX fires CO 2 ) have a weak relationship with the respiration-weighted temperature anomalies (Fig. 8), and a maximum slope of ∼ 0.05 ppm K −1 .

Fossil fuel
The fossil fuel contribution (CTX ff CO 2 ) to the CTX CO 2 contains significant interannual variability (∼ 23 % in Lamont, Bremen and Białystok, and ∼ 10 % in Park Falls, see Table 5). Because we detrend the TCCON data and CT2011 assimilation output of total X CO 2 by the annual measured CO 2 growth rate, some of the variability in the annual fossil fuel emissions will be removed from the X CO 2 anomalies and will not contribute to the X CO 2 −δT slope. Any remaining variability attributed to the fossil fuel signature is likely due partially to the dynamical effects described in Sect. 7.1, which would act to reduce the overall X CO 2 −δT slopes, because fossil fuel emissions primarily occur in the northern midlatitudes, and are of opposite sign to the biospheric drawdown. If we assume that ∼2 PgC are emitted from fossil fuels in the Northern Hemisphere growing season, and that the growing season NEE represents a ∼5 Pg C sink, the fossil fuel emissions represent about one-third of the overall flux. If the dynamics produce about half of the temperature sensitivity, as we suggest in Sect. 7.1, then about one-sixth of the overall signal will be caused by the interaction of fossil fuel emissions and dynamics. This would produce a temperature sensitivity of ∼0.2 ppm K −1 on a hemispheric scale, which is significantly smaller than our uncertainties. The CT2011 fossil fuel signal (CTX ff CO 2 ) is not significantly correlated with the respiration-weighted surface temperature (Fig. 8, bottom  panel), but is not inconsistent with a small, negative slope.
The impact of fossil fuel emissions on TCCON data may be largest in western Europe. Bremen has a much smaller X CO 2 − δT relationship than most of the Northern Hemisphere in both the observations and models. This is a persistent feature of western Europe, which has consistently lower slopes than the rest of the northern midlatitudes (Fig. 7). Further, the rough difference between the Park Falls and Bremen slopes for both the TCCON results and the SiB2009 simulations is ∼1.3 ppm K −1 , suggesting that this difference may reflect the large net efflux from western Europe.

NEE
The most significant contribution to interannual variability in CTX CO 2 is the terrestrial biosphere component, which accounts for 60 % in Bremen, increasing to 82 % in Park Falls (see Table 5). Respiration is directly influenced by surface temperature, and GPP is indirectly influenced by temperature through soil moisture. In SiB, both ecosystem respiration and GPP cause interannual variability in SiBX CO 2 . For example, the 2006 and 2009 SiBX CO 2 seasonal cycle minima are similarly deep compared with 2010, and possess a similar growing season NEE. In 2006, the NEE was more negative relative to 2010 due to a decrease in respiration throughout the growing season. In contrast, in 2009, this was due to an increase in GPP relative to 2010. However, only the interannual variability in the integrated growing season ecosystem respiration is Atmos. Chem. Phys., 13, 9447-9459 significantly correlated with the surface temperature anomalies (R 2 ∼ 0.9, Fig. 9). The interannual variability in growing season GPP is not correlated with δT , suggesting that R e is by far the stronger driver of this X CO 2 −δT relationship in SiB. This is consistent with the results of Vukićević et al. (2003), who show that the temperature sensitivity of GPP is less than that of respiration in all versions of their model. However, "greenness" indices of plant phenology (such as NDVI) have been shown to be poor predictors of GPP in the boreal coniferous region in winter and during cloudy periods (Wang et al., 2004).
NEE depends on soil moisture through GPP and/or R e (Raich et al., 1991;Denning et al., 1996;Zhao and Running, 2010), and so we evaluate soil moisture using the global gridded Palmer Drought Severity Index (PDSI, Dai et al., 2004;Dai, 2011a, b), which is available through 2010. PDSI values that are negative indicate dry conditions, and positive indices indicate wet conditions. Annual averages of PDSI from 30 • N to 60 • N weakly correlate with δT (R 2 = 0.33) for the years 2004-2010. However, it is difficult to determine whether the X CO 2 drawdown values from the TCCON measurements correlate significantly with PDSI. Only the Park Falls dataset has sufficient overlap with the available PDSI values, and it has a negative slope (drier conditions yield higher X CO 2 values). The other sites, even if they possess more than two years of coincident measurements, do not encompass a sufficiently large PDSI range to compute significant slopes.

Summary and future study
Interannual variability in the seasonal cycle minima of column-averaged dry-air mole fractions of CO 2 is correlated with summertime surface temperature anomalies in boreal regions. The CarbonTracker 2011 assimilation and GEOS-Chem simulations suggest that this relationship is caused both by dynamical mixing and biospheric activity, in roughly equal proportion. The effects of interannual variability in emissions from fossil fuel combustion and fires appear to be small and uncorrelated with surface temperature. However, dynamical effects on fossil fuel emissions should have an opposing effect on the X CO 2 −δT relationship. The CT2011 and GEOS-Chem/SiB X CO 2 −δT relationships are generally weaker than observed.
It is clear that there are several avenues worth pursuing to investigate further the processes responsible for the X CO 2 −δT relationship. It is important to try other realizations of the dynamics: using different transport models such as TM5, or different underlying wind fields, such as ECMWF or NCEP. Because our results are derived from columnaveraged atmospheric CO 2 measurements, which are relatively insensitive to planetary boundary layer (local) dynamics, we anticipate that differences in the large-scale dynamics will dominate the variability. A future investigation of the dynamical effects on fossil fuel emission-driven interannual variability should include a GEOS-Chem/SiB run with and without fossil fuel emissions in the model. We anticipate that, in the absence of fossil fuel emissions, the Bremen slopes would be similar to those measured elsewhere. The role of the oceans should be explored in more detail, including an assessment of the impacts of temperature variability and sea ice on gas exchange, and changes in ocean biogeochemistry associated with different climate modes.
To disentangle the effects of NEE on the observed variability further, both the effects of GPP and R e should be probed. Although interannual variability in GPP is not correlated with δT in SiB, Guerlet et al. (2013) have attributed the interannual variability in the GOSAT 2009-2010 seasonal cycles to reduced carbon uptake (GPP) through their inversions. However, the 2010 Northern Hemisphere had exceptional warming, causing record-breaking heat waves throughout eastern Europe and Russia, and fires in western Russia (Barriopedro et al., 2011). The correlation between δT and GPP should be explored in other biospheric models, and by using chlorophyll fluorescence (e.g., Frankenberg et al., 2011) or δ 18 O (e.g., Welp et al., 2011) instead of, or in conjunction with, the current proxies for plant phenology (e.g., NDVI). Further, the sensitivity of R e to temperature (Q 10 ) and moisture should be studied further, through the use of δ 13 C in situ measurements (e.g., Bowling et al., 2002) or other means.

D. Wunch et al.: Covariation of X CO 2 with boreal temperature
Developing robust metrics for respiration and GPP responses to temperature is critical for reducing uncertainties in Earth system models and for diagnosing the vulnerability of permafrost carbon pools to future change.