Oxygen isotope anomaly in tropospheric CO2 and implications for CO2 residence time in the atmosphere and gross primary productivity

The abundance variations of near surface atmospheric CO2 isotopologues (primarily 16O12C16O, 16O13C16O, 17O12C16O, and 18O12C16O) represent an integrated signal from anthropogenic/biogeochemical processes, including fossil fuel burning, biospheric photosynthesis and respiration, hydrospheric isotope exchange with water, and stratospheric photochemistry. Oxygen isotopes, in particular, are affected by the carbon and water cycles. Being a useful tracer that directly probes governing processes in CO2 biogeochemical cycles, Δ17O (=ln(1 + δ17O) − 0.516 × ln(1 + δ18O)) provides an alternative constraint on the strengths of the associated cycles involving CO2. Here, we analyze Δ17O data from four places (Taipei, Taiwan; South China Sea; La Jolla, United States; Jerusalem, Israel) in the northern hemisphere (with a total of 455 measurements) and find a rather narrow range (0.326 ± 0.005‰). A conservative estimate places a lower limit of 345 ± 70 PgC year−1 on the cycling flux between the terrestrial biosphere and atmosphere and infers a residence time of CO2 of 1.9 ± 0.3 years (upper limit) in the atmosphere. A Monte Carlo simulation that takes various plant uptake scenarios into account yields a terrestrial gross primary productivity of 120 ± 30 PgC year−1 and soil invasion of 110 ± 30 PgC year−1, providing a quantitative assessment utilizing the oxygen isotope anomaly for quantifying CO2 cycling.

The net growth rate and level of CO 2 in the atmosphere represents a dynamic balance between anthropogenic activities and natural sources and sinks 1 . The diurnal and seasonal cycles, however, are largely affected by terrestrial photosynthesis and respiration [2][3][4] . The oxygen isotopic composition of atmospheric CO 2 has been shown to be a powerful tracer for improving understanding of carbon and water cycles involving CO 2 (ref. [2][3][4][5][6][7][8][9], providing a unique way to estimate terrestrial gross primary productivity 2,5,8,9 . Oxygen has three stable isotopes ( 16 O, 17 O, and 18 O). The 18 O/ 16 O isotope ratio is used widely for studying aspects of the carbon and water cycles of natural systems 2,3,5-9 , but 17 O/ 16 O has rarely been used owing to added analytical difficulty. Here, we compare and analyze the results of the triple-oxygen isotope composition of surface air CO 2 from northern hemisphere sites in the western Pacific (South China Sea and Taipei, Taiwan) 10 and available data from the middle east (Jerusalem, Israel) 11 and western North America (La Jolla, United States) 12 , to provide deeper insight into the global CO 2 cycle of the atmosphere, biosphere, and hydrosphere (The data from Göttingen, Germany 13 are not included because of the presence of unknown drift in their two-year 17 O data, though their first year of data agree well with the data reported in this work. Our approach is discussed in detail in a following section). We derive the residence time of CO 2 in the atmosphere and gross primary production (GPP) from the integrated data set and discuss how interhemispheric transport affects these quantities. Given that the tropospheric mixing time in each hemisphere is much shorter than the interhemispheric mixing time 14,15 and the latter shorter than the CO 2 residence time derived here (see below), the compiled data should be a valid approximation of the global average.
We show that the terrestrial flux (the CO 2 cycling flux between the atmosphere and terrestrial biosphere) can be quantified using the Δ 17 O values of CO 2 , where

Sources and processes defining Δ 17 O values
The size of the anomaly is dependent upon the choice of λ and the reference scale. Here, we choose 0.516 for the slope and the most commonly used scale for oxygen, V-SMOW. We note that the selection of reference scale does not affect interpretation, provided the variation and partitioning among the three oxygen isotopes are properly accounted for. Given that the carbon flux estimation presented in this paper is based on the deviations of the oxygen anomalies of reservoirs/processes from that measured in atmospheric CO 2 , it is natural to take the λ value describing the variation of the triple-oxygen isotopic partitioning in tropospheric CO 2 . Processes that affect CO 2 isotopologues in the troposphere are terrestrial, oceanic and anthropogenic, with the first being dominant.
In the terrestrial biosphere, leaf water transpiration governs the variation of oxygen isotopes in CO 2 ; the mean, however, is largely determined by water-CO 2 equilibration catalyzed by carbonic anhydrase 2,3,[5][6][7]9,33 . It has been found previously that the transpiration λ value is a function of air relative humidity 21 , whereas dependence on other meteorological variables such as temperature and soil water isotopic composition has not been observed. We set λ = 0.516, as it represents the transpiration λ at 75% relative humidity, a globally averaged humidity near the surface 34 . This λ value is essentially the same as that of CO 2 we obtained (0.518 ± 0.004, see below) for the western Pacific, which had an average relative humidity of 76 ± 4% during 2010-2015 (data obtained from the Center for Weather Bureau, Taiwan; site code: 466920; the value changes slightly to 72 ± 11% if considering only day time between 6 AM to 6 PM). However, given the sparse spatiotemporal coverage of the data, the governing slope cannot be firmly decided. More data taken under a variety of environmental conditions are needed to set a better constraint. From the current understanding of the processes occurring, we consider plant transpiration the most important process affecting the variation of atmopsheric CO 2 (e.g., see Landais et al. 21 and Cuntz et al. 3 and references contained therein). We stress that the value of λ does not affect the flux interpretation shown below, as long as equation (2) is used, but the selection must best represent the variation of atmospheric CO 2 .
A schematic diagram (not to scale) that describes various sources and processes modifying CO 2 isotopologues is shown in Fig. 1, which summarizes the oxygen isotope transport at steady state. On a global scale, equilibrium processes are the major controllers in oxygen isotope dynamics; we show below in the Box model section that kinetic fractionations are insignificant. (Previous work 13 utilized a time-dependent model showing a measurable seasonal cycle of Δ 17 O at an amplitude of ~0.05‰. The present work assesses the global carbon budget at steady state, leaving the assessment of spatiotemporal variability, including seasonality, to a latter paper when data covering a variety of spatial and temporal scales are available). Three sources/processes are considered: terrestrial (meteoric) water, ocean water, and anthropogenic CO 2 , with the last inheriting the atmospheric O 2 isotopic composition. The variation of the triple oxygen isotopic composition of meteoric water follows a slope of 0.528 (ref. 20 ). Subsequent isotopic exchange between water and atmospheric CO 2 modifies the oxygen isotopic composition of CO 2 following a slope of 0.5229 (ref. 11 ). In addition, in the terrestrial biosphere, plant transpiration changes the source water following a slope of 0.516, a value chosen to represent the average slope at the globally averaged relative humidity of 75% (ref. 21,34 ), as discussed above. These three values (0.528, 0.5229, and 0.516) determine the Δ 17 O value of CO 2 mediated by processes involving water. We use the scheme in Fig. 1 to follow changes in oxygen isotopic composition, starting with meteoric water (point A). Transpiration that affects leaf water changes source meteoric water from A to B; the line AB follows the slope of 0.516. Subsequent water-CO 2 isotopic exchange determines the composition of oxygen in leaf CO 2 , following the slope of 0.5229 to point C. For exchange with ocean water, only water-CO 2 equilibration is involved, that changes CO 2 to the blue "X. " When the CO 2 enters the stratosphere, the coupled CO 2 -O 2 -O 3 photochemistry moves the tropospheric CO 2 to D; this path goes through mean tropospheric CO 2 (diamond) and has slope of >1 (ref. [29][30][31][32] ). Large-scale circulation and synoptic eddy mixing bring the modified CO 2 back to the troposphere. The tropospheric CO 2 (diamond) represents a balance among stratospheric (grey), terrestrial (yellow), oceanic (blue "X"), and the final component -anthropogenic CO 2 (red symbol). Figure 1B, plotting Δ 17 O versus ln(1 + δ 18 O), shows that our reference λ is chosen (0.516) so that plant transpiration does not change the Δ 17 O value of terrestrial water. This is a rotation of Fig. 1A such that a slope of 0.516 becomes 0, and this removes the dominant variation along the correlation in Fig  In addition to well-quantified photochemical processes in the middle atmosphere, there are two known processes that modify Δ 17 O: combustion and isotope exchange with water. The former produces CO 2 with Δ 17 O = −0.21‰, the same as air O 2 (ref. 19 ), used for the anthropogenic component below.
Isotopic exchange with water can be estimated using water-CO 2 equilibrium 11,23 . Under the assumed definition of λ derived from equation (1), Δ 17 O of CO 2 emitted from sources involving water-CO 2 equilibration processes like respiration and soil invasion, following the slope of 0.5229 (ref. 11 34 ). Globally averaged meteoric water has Δ 17 O mw = −0.046 ± 0.005‰ (1 standard error; or 0.032 ± 0.003‰ at λ = 0.528; ref. 20 ), excluding highly depleted waters having δ 18 O less than −10‰ in high latitude regions covered by snow and/or ice. (Here, standard error represents the error of a sample mean; standard deviation describes the error of a single measurement, the spread of replicate analyses of a single sample, or the spread of an ensemble). Figure 3 shows that the values of Δ 17 O mw from various regions are comparable.
The other largest water reservoir is the oceans. The globally averaged 20   is 0.284 ± 0.001‰. The observation that Δ 17 O of meteoric water is lower than that of oceanic water suggests that the Δ 17 O approach has greater sensitivity to terrestrial processes than to oceanic ones.
Below we combine the values calculated above with data from several locations around the world (La Jolla, CA, USA; Jerusalem, Israel; Taipei, Taiwan; and the South China Sea) to put constraints on GPP and the oxygen isotope residence time of CO 2 in the atmosphere. The oxygen isotope residence time is defined by the ratio of the atmospheric CO 2 mass loading (M) and the CO 2 mass flux between the atmosphere and biosphere/hydrosphere (F sur ). The flux is inferred from the mass balance calculation obtained using the triple oxygen isotopic composition of tropospheric CO 2 shown below.

Methods
In addition to using data available in the literature from Jerusalem 11 , La Jolla 12 , and western Pacific regions 10 , we have continuously collected air for isotopic analysis of CO 2 in three locations: (1) Academia Sinica campus (abbreviated AS; 121°36′51″ E, 25°02′27″ N; ~10 m above ground level or 60 m above sea level) in Taipei, Taiwan and (2) the campus of National Taiwan University (NTU; 121°32′21″ E, 25°00′53″ N; ~10 m above ground level or 20 m above sea level; ~10 km southwest of Academia Sinica). To check the reported Δ 17 O values in the eastern Pacific 12 , we have also collected and analyzed CO 2 from Los Angeles, California at a latitude slightly higher than La Jolla, along the coast on Palos Verdes peninsula (118°10.9′ W, 33° 44.7′ N; PVD). Data reported in this work and analyzed in Taiwan are provided in full in the supplementary material.
Analytical methods are described in detail elsewhere 10,18,[36][37][38] and summarized here. Air from western Pacific regions for isotope analysis was collected in pre-conditioned 1-liter Pyrex bottles, achieved by passing dry, high purity nitrogen through the bottles overnight. The sampling bottles used for concentration (~350-ml bottle) and isotope (1-liter) analyses were connected in series. Samples were collected and compressed to 2-bar after flushing the bottles for 5 minutes with ambient air at a flow rate of ~2 liter per min. Moisture was removed during sampling using magnesium perchlorate to minimize subsequent isotope exchange between CO 2 and water 39 . Concentration of CO 2 is measured using a LI-COR infrared gas analyzer (model 840A, LI-COR, USA), with reproducibility better than 1 ppmv. The PVD samples were collected on Saturday afternoons at about 14:00 PST, into 2-liter evacuated Pyrex flasks after passing through Mg(ClO 4 ) 2 . Carbon dioxide was separated from the air samples cryogenically and measured, following the method described in Newman et al. 38 . The CO 2 -O 2 oxygen isotope exchange method developed previously 36,37 was used to measure the Δ 17 O of CO 2 samples. Isotopic analyses were done using a FINNIGAN MAT 253 mass spectrometer in the dual inlet mode. The analytical precision obtained for a single measurement of the Δ 17 O value of CO 2 is better than 0.01‰ (1-σ standard deviation).

Results
The concentrations of the isotopologues of CO 2 at South China Sea (SCS) are close to those reported at Mauna Loa measured by the National Oceanic and Atmospheric Administration -Earth System Research Laboratory (data available online at http://www.esrl.noaa.gov/gmd/dv/iadv/), with [CO 2 ] 395.4 ± 7.3 ppmv and δ 13 C −8.47 ± 0.22‰ (1-σ standard deviation). In Taipei (AS + NTU), the averaged [CO 2 ] is 416.2 ± 18.3 ppmv, δ 13 C is −9.22 ± 0.83‰, and δ 18 O is 40.65 ± 0.82‰. [CO 2 ] varies between ~350 and 475 ppmv, with low values during the day and high at night, representing the combined effect of natural biogeochemical cycle (photosynthesis and respiration), anthropogenic emissions, and boundary mixed layer height diurnal variation 38 . During day time hours, photosynthesis dominates, resulting in reduction of CO 2 concentration and less negative δ 13 C and δ 18 O values. The CO 2 content is lower during the day than night also due to dilution as the boundary layer deepens.
The mean values obtained for [CO 2 ] and δ 13 C, as compared to SCS, show that in Taipei a clear contribution from anthropogenic emissions is seen. Given the proximity of the AS and NTU stations (~10 km apart), air transport time is shorter than the CO 2 oxygen isotope residence time shown below, resulting in, on a yearly basis, similar levels of CO 2 isotopologues, including Δ 17 O. Figure 2A compiles tropospheric CO 2 data in a plot of ln(1 + δ 17 O) vs. ln(1 + δ 18 O). The least square linear regression analysis of the data obtained in Taipei yields a slope of 0.518 ± 0.004 (excluding some outliers at δ 18 O < 38.5‰) and intercept of 0.3 ± 0.2‰ (1 standard error, R 2 = 0.99). The observation that the same slope is obtained as for transpiration at ~75% relative humidity (close to the value at the sampling sites in western Pacific, 76 ± 4% averaged between 2010 and 2015) suggests that transpiration is likely a controlling process affecting the variation of the triple oxygen isotopic composition of near surface CO 2 . As the variation of oxygen isotopic compositions in CO 2 at the two Taipei sampling sites is biogeochemically mediated, one may use Δ 17 O to estimate the actual "flux" between the atmosphere and soil/leaf, which in turn gives the value for GPP. The overall ratio of ln(1 + δ 17 O)/ln(1 + δ 18 O) for Taipei is 0.524 ± 0.001, consistent with that of a water-CO 2 equilibrium value of 0.5229 ± 0.0001 (ref. 11,23 ), further verifying that the oxygen of near-surface CO 2 in this region is primarily affected by biogeochemistry with minor influences from the stratosphere and human activities. The Δ 17 O values in Taipei vary from 0.216 to 0.415‰ (Fig. 2B) with an average of 0.335 ± 0.039‰ (1-σ standard deviation of the range), a value similar to 0.31 ± 0.06‰ obtained in La Jolla 12 and 0.321 ± 0.007‰ in Jerusalem 11 . (Note that the values for the last two have been re-scaled with respect to λ = 0.516, for the sake of consistency among the data sets). The value at SCS is 0.335 ± 0.033‰. The value reported at La Jolla is the lowest among the four. We then check the possibility that the difference between La Jolla and Taiwan is caused by the difference in Δ 17 O scale between the two labs. By taking CO 2 samples collected also in the eastern Pacific region but at a higher latitude near a southwest facing beach on the Palos Verdes peninsula (PVD). The obtained values for 2015 are 0.317 ± 0.032‰, essentially the same as the decadal mean from La Jolla. The average of the mean Δ 17 O values from La Jolla, Jerusalem, Taipei, and the South China Sea (0.326 ± 0.005‰) is used below in the box model calculation of GPP and oxygen isotope residence time. The averaged Δ 17 O value decreases by only 0.004‰ if just the two largest datasets, from Taiwan and La Jolla, are averaged. Therefore, we note that whether the SCS, Jerusalem, and PVD data are included does not change the conclusion presented below. For example, a 0.01‰ reduction in Δ 17 O results in ~0.2-year decrease in the oxygen isotope residence time. Moreover, given that the oxygen isotope residence time is on the order of one year only and sources and processes that change the isotopologues of the tropospheric CO 2 are variable spatially and temporally, spatiotemporal inhomogeneity in Δ 17 O is expected. The average of the mean values from more locations should be more representative for the global Δ 17 O.

Box model
A box model is employed to assess various contributing processes for Δ 17 O. At steady state, the mass balance equation for δ (where δ is either δ 17 O or δ 18 O), following Cuntz et al. 3 , can be written as follows: where C is the volume mixing ratio of CO 2 , M is the mass of the atmosphere, subscripts "a", "st", "anth", "l", "r", "s", and "o" of δ's represent the δ values of the sampled air, the stratosphere, anthropogenic emissions, and leaf, respiration, soil, and ocean water, respectively, ε's are the associated kinetic fractionation factors, and F is the flux in and out of a reservoir such that the subscript "la" refers to leaf-to-air, "al" air-to-leaf, and "ao" air-to-ocean. We then rewrite the equation (3) in terms of Δ in steady state as follows. where the kinetic terms ((F la − F al ) × ε l × (λ l − λ 0 ) + F r × ε s × (λ s − λ 0 ) + F ao × ε o × (λ o − λ 0 )) become negligible (even with extreme values for F at 500 PgC year −1 or λ l,s,o at 0.529, the value for equilibrium water between condensed and vapor phases, the isoflux is found to be less than 1‰ PgC year −1 ). We note that Δ in the linear definition in equation (4) (4) results in an error about 10% in each term derived, i.e., ~40 PgC year −1 biased too high in F la + F r + F s and ~0.2 year too short in the resident time of CO 2 in the atmosphere (though still within the error of the estimation). Parameters and values used in the box modeling are summarized in Table 1. . We further take F la ≈ 2F r ≈ 2F s (from the fact that global net productivity is much less than the gross productivity, and the assumption of catalyzed soil invasion 7 ), and with this, the terrestrial flux F land equal to F la + F r + F s . Then equation (4) can be reduced to  16,42 ). Figure 4 summarizes the derived terrestrial CO 2 flux and residence time in the atmosphere; in this particular model, the cross-hemispheric transport and mixing are not included, as the hemispheric difference in Δ 17 O was predicted to be small (<0.01‰; ref. 16 ). Sensitivities of the derived quantities with respect to the variations of the relative importance of ocean flux, cross-tropopause exchange flux, soil invasion, and Δ 17 O value in the southern hemisphere are presented in Fig. 5. See below for the detail on the error assessment.
The current steady state box model is an updated version of Hoag et al. 16 ; the major surface resetting processes are included explicitly to distinguish the terrestrial (re-)cycling fluxes from the oceanic. Previous works (Cuntz et al. 3 and Hofmann et al. 13 ) solve time-dependent equations (equation 3). In this case, kinetic fractionations becomes important, and spatiotemporal variabilities in all components of carbon/oxygen cycling and recycling fluxes are significant. As a result, spatial and temporal inhomogeneity in sampling has to be considered and evaluated critically. When more data are available, natural variability in carbon fluxes can be assessed. Therefore, for examining the global carbon flux at steady state, we take measurements spanning as much space and time as possible. For error analysis, a standard error from each component of measurements is adopted, leaving standard deviation for representing spatiotemporal variability assessment in time-dependent models.

Assessing gross primary productivity
Plant uptake scenarios affect the estimates of GPP and soil invasion. GPP can be estimated as follows: al la eq c land s where θ eq represents the degree of hydration of CO 2 in stomata and κ c is a measure of stomatal conductance which can be expressed by where C c is the CO 2 concentration in chloroplasts at the site of CO 2 hydration and C a is the atmospheric concentration. The factor 0.88 is used to account for leaf respiration 43 . For C 3 plants, C c /C a = 2/3; for C 4 plants, C c /C a = 1/3, assuming that C c is equal to intracellular CO 2 concentration 44 . A generally averaged C c /C a is 0.57 or κ c = 1.33 (ref. 8,9 ). Overall, the value of κ c varies between 1.33 and 2.93 (ref. 3,9 ). θ eq is also variable 33,45 . Careful error assessment is made to quantify the effect of this variation on the value of GPP (Fig. 5C). See Table 1 for the summary of the parameters and the associated errors and variations used in this work. Errors resulted from incomplete understanding of sources and the Δ 17 O values in various water bodies (meteoric, soil, and oceanic waters) and global atmospheric CO 2 (i.e., tropospheric and stratosphere-to-troposphere stream Δ 17 O). The corresponding errors are summarized in Table 1. There are, however, other related processes whose λ values remain poorly known, of which transpiration and plant uptake scenarios (i.e., chloroplast CO 2 concentration and enzyme catalyzed oxygen exchange efficiency) are the most uncertain. Transpiration is sensitive to plant species and ambient air relative humidity 21 . If we assume an average considered in Landais et al. 21 , the variation in λ results in a variation of 0.017‰ for Δ 17 O for a 0.05 change in relative humidity; this term contributes ~30% of the error (~1/3 from atmospheric CO 2 measurements and the remaining from the cross-tropopause exchange Δ 17 O flux) in estimating the global oxygen isotope residence time and terrestrial carbon flux. Proper error propagation is made to assess the error in the terrestrial flux (F land ) derived in equation (5), followed by a Monte Carlo approach to evaluate the errors from κ c and θ eq in equation (6) for GPP estimate.
A function f with n independent variables is expressed as:  The error of each quantity x i is given by the standard deviation (σ i ) or standard error σ ( ) i . The former refers to the error of a single measurement while the latter is used to represent the error of a sample mean. With a large set of measurements, we use σ i to represent the error of the sample mean x i . The errors in equation (4) refer to such an error. Following standard error propagation, the error of the function f is , if we measure them once, the error, or sample inhomogeneity due to inhomogeneity in sources and processes, increases to 600 PgC year −1 . For errors resulting from global carbon cycling scenarios, a Monte Carlo simulation is employed. The final error is estimated from 500 × 500 × 500 = 125 M (500 samplings for the terrestrial flux and 500 each for two variables in the carbon model in equation (6); the errors are converged with >300 random samplings) calculations of the model using ) is a factor of 10 higher than the oceanic isoflux (5‰ PgC year −1 ), suggesting that the Δ 17 O approach has greater sensitivity to terrestrial processes than oceanic processes. For easier comparison, we have multiplied the oceanic and terrestrial isofluxes by −1.
repeated random sampling of each of the three variables. From the literature, κ c ranges from 1.33 to 2.97 (depending critically on C c ; ref. 3,9 ) and θ eq from 0.7 to 0.78 (affected by enzyme carbonic anhydrase activity in C 4 plants; ref. 9,33,44 ). To consider final total errors in GPP and F s , we perform Monte Carlo simulations to randomly determine values for κ c and θ eq over the aforementioned ranges. The distribution of these two parameters is assumed to be uniform, with the measurement error of Δ 17 O assumed to be following normal distribution. Figure 4 summarizes the model results calculated using equation (5) and a value for Δ 17 O of 0.326 ± 0.005‰ for the troposphere, the average for the four locations discussed above. The current mass loading of atmospheric CO 2 (M) is 828 ± 10 PgC (ref. 1,45 ). The CO 2 oxygen isotope residence time τ is given by M/F sur , where F sur = F land + F ao . Taking into account the aforementioned uncertainties of the parameters in equation (5), the terrestrial flux F land is determined to be 345 ± 70 PgC year −1 , and τ is 1.9 ± 0.3 years, consistent with previous estimates 2, 3,8,9 . The estimate is insensitive to the partitioning between ocean and terrestrial fluxes because of the sensitivity of the Δ 17 O approach to the terrestrial processes (Fig. 4). However, we show below that we cannot constrain the value for GPP better than other methods, because of unknown quantities for soil invasion and degree of isotopic equilibrium between leaf water and stomatal CO 2 .

Oxygen isotope residence time and gross primary productivity
No significant advancement toward quantifying soil invasion has been made since Wingate et al. 7 . The reported flux can be as low as <10 PgC year −1 (ref. 6 ) to as high as 450 PgC year −1 (ref. 7 ), depending on the catalyzed hydration activity (via enzyme carbonic anhydrase). Recently, the hydration activity has been found to likely be high, resulting in an invasion flux as high as respiration 7 , and we choose this as our best estimate, i.e., F s = F r . By definition, GPP is the sum of respiration F r and NEP (net ecosystem productivity). So equation (6) can be expressed by The globally estimated NEP is 10 PgC year −1 (ref. 46 ). Once the values of κ c and θ eq are chosen, along with the value of F land reported above, GPP and F s can be derived.
With a previously suggested plant uptake scenario (κ c = 1.33 and θ eq = 0.78; ref. 8 ) and an independent constraint for net ecosystem productivity 46,47 , we derive an estimate of 130 ± 25 PgC year −1 for GPP, with a best guess for soil invasion of 120 ± 20 PgC year −1 (calculated from equation (10)). A Monte Carlo simulation that considers various carbon cycling models 3,8,9,33,44 , including plant types and degree of oxygen equilibrium with various water bodies in the biosphere and hydrosphere, gives the estimates of GPP and soil invasion to 120 ± 30 and 110 ± 30 PgC year −1 , respectively. The estimated GPP is toward the lower end of Welp et al. 8 but close to that of Beer et al. 4 . Given that the value of Δ 17 O is sensitive to the isofluxes between atmospheric CO 2 and water bodies, we expect that extended studies with multiple CO 2 isotopologues into C 3 -dominated regions such as the Amazonian rainforest, C 4 -dominated lands such as grasslands in North America, and vegetation-sparse areas such as the Canadian arctic could have great potential to refine the partitioning between photosynthesis, respiration and soil invasion, and thus to provide a better estimate of the terrestrial GPP.  Figure 5D shows that the GPP and residence time are sensitive to the Δ 17 O value in the southern hemisphere; the sensitivity is obtained by assuming F s = 110 PgC year −1 , F land = 345 PgC year −1 , κ c = 1.33, and θ eq = 0.78. If Δ 17 O in the southern hemisphere is 0.01‰ (0.02‰) higher than that in the northern hemisphere, the GPP value increases to 170 (200) PgC year −1 ; the residence time then reduces to 1.7 (1.5) years. As soon as the value of Δ 17 O in the southern hemisphere is measured, the global residence time can be determined and GPP can be better quantified.
In summary, the triple-oxygen isotopic composition of CO 2 constrains the global oxygen isotopic residence time in the atmosphere to 1.9 ± 0.3 years, compared to 0.9-1.7 years (ref. 8 ) or longer 2,3,9 reported previously. The terrestrial flux is quantified to be 345 ± 70 PgC year −1 , falling in the range reported in the literature, 200-660 PgC year −1 (ref. 2,3,8,9 ). Because of the isotope recycling time of CO 2 , the spatial inhomogeneity of Δ 17 O obtained between localities shows that the commonly used δ values can be applied to Δ 17 O to refine knowledge of the flux partitioned between respiration/soil invasion, photosynthesis, and air-sea exchange. CO 2 sampling campaigns in the remote Pacific and southern hemisphere oceans can better remove interference from terrestrial processes, to quantify the oceanic flux. High-resolution global and regional models assimilating CO 2 isotopologues (Δ 17 O in particular) with online carbon and water cycle modules can potentially strengthen our understanding of the associated processes at molecular scales. We expect that existing models 3 coupled with a cross-tropopause exchange module extending surface biogeochemical models to include stratospheric processes will greatly improve our estimates and provide extraordinary precision to probe the associated fluxes in the global carbon and water cycles involving CO 2 .