Simulation of carbon isotope discrimination of the terrestrial biosphere

We introduce a multistage model of carbon isotope discrimination during C3 photosynthesis and global maps of C3/C4 plant ratios to an ecophysiological model of the terrestrial biosphere (SiB2) in order to predict the carbon isotope ratios of terrestrial plant carbon globally at a 1° resolution. The model is driven by observed meteorology from the European Centre for Medium‐Range Weather Forecasts (ECMWF), constrained by satellite‐derived Normalized Difference Vegetation Index (NDVI) and run for the years 1983–1993. Modeled mean annual C3 discrimination during this period is 19.2‰; total mean annual discrimination by the terrestrial biosphere (C3 and C4 plants) is 15.9‰. We test simulation results in three ways. First, we compare the modeled response of C3 discrimination to changes in physiological stress, including daily variations in vapor pressure deficit (vpd) and monthly variations in precipitation, to observed changes in discrimination inferred from Keeling plot intercepts. Second, we compare mean δ13C ratios from selected biomes (Broadleaf, Temperate Broadleaf, Temperate Conifer, and Boreal) to the observed values from Keeling plots at these biomes. Third, we compare simulated zonal δ13C ratios in the Northern Hemisphere (20°N to 60°N) to values predicted from high‐frequency variations in measured atmospheric CO2 and δ13C from terrestrially dominated sites within the NOAA‐Globalview flask network. The modeled response to changes in vapor pressure deficit compares favorably to observations. Simulated discrimination in tropical forests of the Amazon basin is less sensitive to changes in monthly precipitation than is suggested by some observations. Mean model δ13C ratios for Broadleaf, Temperate Broadleaf, Temperate Conifer, and Boreal biomes compare well with the few measurements available; however, there is more variability in observations than in the simulation, and modeled δ13C values for tropical forests are heavy relative to observations. Simulated zonal δ13C ratios in the Northern Hemisphere capture patterns of zonal δ13C inferred from atmospheric measurements better than previous investigations. Finally, there is still a need for additional constraints to verify that carbon isotope models behave as expected.


Introduction
[2] Atmospheric CO 2 concentrations have been increasing for the past 150 years due to combustion of fossil fuels, manufacture of cement, and biomass burning. In the 1980s, this anthropogenic CO 2 was added to the atmosphere at a rate of approximately 5.4 ± 0.3 Pg-C/year; in the 1990s the rate was 6.3 ± 0.4 Pg-C/year [Intergovernmental Panel on Climate Change (IPCC), 2001]. However, over these same periods, atmospheric CO 2 only increased at rates of [3] Carbon isotopes can be helpful in investigations of the sink because they give information about both sources and processes of CO 2 exchange. At the local level, carbon isotopes are helpful in understanding (1) plant water use efficiency and the response of plants to changes in precipitation and relative humidity [Farquhar and Richards, 1984;Farquhar et al., 1988;Condon et al., 1993;Hall et al., 1993;Schulze et al., 1996Schulze et al., , 1998Ekblad and Högberg, 2001;Bowling et al., 2002;Ometto et al., 2002], (2) variation in light distribution and stand structure [Berry et al., 1997;Buchmann et al., 1997aBuchmann et al., , 1997bLe Roux et al., 2001], (3) recycling of respired CO 2 [Keeling, 1961;Schleser and Jayasekera, 1985;Sternberg, 1989;Lloyd et al., 1996;Sternberg et al., 1997], and (4) determining the relative contributions of photosynthesis and respiration to total net ecosystem exchange [Yakir and Wang, 1996;Bowling et al., 2001]. At the global scale, carbon isotopes have been used to (1) determine spatial and temporal distribution of sources and sinks and, in some cases, determine partitioning of oceanic and terrestrial fluxes [Francey, 1985;Keeling et al., 1989Keeling et al., , 1995Tans et al., 1993;Ciais et al., 1995aCiais et al., , 1995bEnting et al., 1995;Fung et al., 1997;Rayner et al., 1999;Bousquet et al., 1999aBousquet et al., , 1999b and (2) to examine relative contributions of C3 and C4 plants to total primary productivity [Lloyd and Farquhar, 1994;Ehleringer et al., 1997;Sage et al., 1999]. Many of these studies show that plant discrimination can change rapidly and that it is very sensitive to environmental changes; however, most of the global models must make simplifying assumptions concerning carbon isotope discrimination of the terrestrial biosphere, and few models function on multiple spatial and temporal scales, so that the results can be directly compared to observations at both flux towers and flask stations.
[4] Previous studies that modeled carbon isotope discrimination of the terrestrial biosphere used monthly mean values of temperature, precipitation, and humidity, were driven by meteorology produced in a GCM, or were produced using modeled vegetation [Lloyd and Farquhar, 1994;Fung et al., 1997;Kaplan et al., 2002], and none operate within the framework of GCM which is simultaneously calculating turbulent exchange of heat, water, and carbon fluxes between the biosphere and atmosphere. We wanted to construct a model that was capable of being driven by observed meteorology and vegetation so that results could be compared to historical records. We also wanted a model that could function online with a GCM so that we could examine systematic relationships between Figure 1. Flow of data in the simulation. ECMWF provides 6-hourly (1) wind speed at 10 m height, (2) temperature at 2 m, (3) relative humidity at 2 m, (4) total incident solar radiation, and (5) precipitation. Monthly vegetation maps (linearly interpolated each day) are from processed satellite data. C3/C4 plant distributions are from ecosystem modeling, satellite data, and maps of agriculture . SiB2 calculates (1) rates of assimilation, (2) leaf boundary layer, stomatal and mesophyll conductance, and (3) CO 2 concentrations in the plant (C s , C i , C c ) and canopy (C a ). C3/C4 distributions are from ecosystem modeling, satellite data, and maps of agriculture ] (see Figure 2 in this paper). CRAX calculates carbon isotope discrimination of C3 and C4 plants. The dashed arrow indicates that CFRAX can be coupled interactively with SiB2 to provide d 13 C of respiration and canopy CO 2 during each time step. That feature is not used in this study. biology and atmospheric dynamics and so that it could be used as a predictive tool. Finally, we wanted a model that would function on local to global spatial scales and hourly to interannual timescales so that spatial and temporal scales of the simulations would be comparable to the spatial and temporal scales represented by data collected at both flux towers and flask stations.
[5] The intent of this paper is to (1) present the mechanics of a model and simulation of terrestrial carbon isotope discrimination, (2) present the results of a global simulation at 1°resolution that is driven by observed meteorology and satellite-derived vegetation characteristics, and (3) evaluate simulation results by comparing them to observations taken at different temporal and spatial scales. Specifically, we compare (1) the simulated relationship between d 13 C of assimilated carbon and mean daily vapor pressure deficit (vpd) to field data from Bowling et al. [2002]; (2) the simulated relationship between d 13 C of assimilated carbon and monthly precipitation to field data from Ometto et al. [2002]; (3) mean d 13 C ratios for four biomes (Broadleaf, Temperate Broadleaf, Temperate Conifer, and Boreal) to observed d 13 C of ecosystem respiration from those biomes [Pataki et al., 2003]; and (4) zonal mean assimilationweighted d 13 C ratios of terrestrial biomass for 20°N to 60°N to d 13 C ratios of CO 2 fluxes inferred for terrestrially dominated stations of the NOAA Globalview Flask Network in the same latitude band [Miller et al., 2003].

Method
[6] SiB2 is a model of the terrestrial biosphere that calculates surface fluxes of sensible and latent heat, radiation, moisture, CO 2 , and momentum for vegetated land points [Sellers et al., 1986Denning et al., 1996a;Schaefer et al., 2002]. SiB2 functions at local, regional, and global spatial scales, and can be driven by observed meteorology or coupled to a global climate model. In the present simulation, SiB2 is driven by assimilated meteorology for 1983 -1993 provided by the European Centre for Medium-Range Weather Forecasting (ECMWF) with a 1°s patial resolution and a 6-hour time step (Figure 1) (see Schaefer et al. [2002] for more results from this particular simulation). ECMWF provides (1) wind speed at 10 m height, (2) temperature at 2 m, (3) relative humidity at 2 m, (4) total incident solar radiation, and (5) precipitation [Gibson et al., 1999]. Time-varying phenological properties such as leaf area index (LAI), the fraction of photosynthetically active radiation absorbed by the green canopy Figure 2. Distribution of C3 and C4 plants  is determined using remote sensing products [Defries and Townshend, 1994;DeFries et al., 1998DeFries et al., , 1999aDeFries et al., , 1999bDeFries et al., , 2000, physiological modeling [Collatz et al., 1998], crop fraction maps [Ramankutty and Foley, 1998], and areal coverage data on crop types [FAO, 1998]. See color version of this figure at back of this issue.
(FPAR), and canopy greenness fraction are linearly interpolated each day from monthly maps of processed normalized difference vegetation index (NDVI) calculated from AVHRR red and near-infrared reflectance [Los, 1998;Los et al., 2000]. Vegetation type is from Defries and Townshend [1994]. Other constant parameters, including canopy height and soil properties, are from look-up tables [Zobler, 1986;Sellers et al., 1996b]. The distribution of C3 and C4 plants ( Figure 2)  is determined using remote sensing products [Defries and Townshend, 1994;DeFries et al., 1998DeFries et al., , 1999aDeFries et al., , 1999b, physiological modeling [Collatz et al., 1998], crop fraction maps [Ramankutty and Foley, 1998], and areal coverage data on crop types [Food and Agriculture Organization (FAO), 1998]. Spatial resolution in the simulation is 1°Â 1°. Recent improvements to SiB2 include the introduction of a six-layer soil temperature submodel based on the work of Bonan [1996Bonan [ , 1998] and a revised surface energy budget that includes prognostic temperature and moisture of the canopy air space reservoir . The calculations are implicit and allow for canopy storage, i.e., memory. Model results have been evaluated at local scales , regional scales Nicholls et al., 2004] and global scales [Denning et al., 1996b].
[7] Another major change to SiB2 is the introduction of a scheme for calculating carbon isotope ratios of terrestrial CO 2 fluxes. Discrimination against 13 C during C3 photosynthesis is the result of a multistage process (Figure 3) involving relatively small isotope effects during transport of CO 2 from the canopy to the chloroplast and a large isotope effect associated with fixation by ribulose bisphosphate carboxylase/oxygenase (rubisco). Similar to equation (A6) in the Appendix of Farquhar and Richards [1984], we model this with three transport steps: (1) molecular diffusion across a laminar boundary layer at the leaf surface, (2) molecular diffusion through a stomatal pore into the stomatal cavity, and (3) dissolution into mesophyll water and aqueous transport to the chloroplast. This is followed by fixation with rubisco or other enzymes. Unlike Farquhar and Richards [1984], we do not include isotope effects associated with either dark respiration or photorespiration. This is because these effects are poorly characterized and probably quite small [Farquhar et al., 1989;O'Leary, 1993;Brugnoli and Farquhar, 2000]. During photosynthesis, a CO 2 concentration gradient from the atmosphere to the chloroplast produces a flux into the plant. At each step, there is a slightly different CO 2 concentration which we designate as follows: C a is CO 2 in the canopy air, C s is CO 2 in the leaf boundary layer, C i is CO 2 in the stomatal cavity, and C c is CO 2 at the chloroplast. A resistance modulates transport from one step to the next: r b for transport from the canopy air into the leaf boundary layer, r s for transport through the stomatal pore, and r c for transport to the chloroplast. These steps are also associated with isotope effects discriminating against 13 C: D b , D s , D diss , and D aq , respectively. Finally, there is an isotope effect during enzymatic fixation with rubisco, D f . The D values that we use are D b = 2.9%, D s = 4.4%, D diss = 1.1%, D aq = 0.7%, Figure 3. Schematic representation of carbon flow through C3 plants in CFRAX. C a , C s , C i , and C c are CO 2 partial pressures in the canopy, at the leaf surface, and in the stomatal cavity and chloroplast, respectively; r b , r s , and r m are resistances to molecular diffusion of CO 2 across the leaf boundary layer, through the stomatal pore and from the stomatal cavity to the chloroplast, respectively; and D b , D s , D diss , D aq , and D f are kinetic isotope effects accompanying molecular diffusion of CO 2 across the leaf boundary layer, through the stomatal pore, and from the stomatal cavity to the chloroplast (dissolution plus aqueous transport) and enzymatic fixation, respectively [Craig, 1953;Mook et al., 1974;Farquhar, 1983;O'Leary, 1984]. and D f = 28.2%, where D is defined as follows. For the reaction A ! B, where R A = [ 13 C/ 12 C] A d A = (R A /R PDB À 1) Â 1000 and R PDB , the 13 C/ 12 C ratio of Pee Dee Belemnite, is equal to 0.0112372 [Craig, 1957]. Units for d and D are per mil (%).
. D s and D b are from theoretical calculations of molecular diffusion of CO 2 through air and across a laminar boundary layer [Craig, 1953;Farquhar, 1983]. D diss and D aq are from laboratory measurements [Mook et al., 1974;O'Leary, 1984]. At the chloroplast, CO 2 can react with rubisco or other enzymes. The fraction reacting with other enzymes is not well known. Raven and Farquhar [1990] have suggested that the upper limit is about 10%. In our simulation, D f is calculated by assuming that 92.5% of the assimilated CO 2 reacts with rubisco with an isotope effect of 30.0% [Brugnoli and Farquhar, 2000], while the remaining reacts with phosphoenyl pyruvate carboxylase (PEPC), with an isotope effect of 5.7% [Farquhar, 1983], giving a net isotope effect (D f ) of 28.2%. Total isotope fractionation during C3 photosynthesis (DPS C3 ) is given by sum of the concentration gradient-weighted isotope effects for each stage, [8] In the present simulation, C a is held constant at 35 Pa, though it can also be allowed to fluctuate naturally in order to isolate the effects of a dynamic canopy air space on d 13 C of biospheric CO 2 fluxes. Net assimilation (A n , where A n is gross photosynthesis minus aboveground autotrophic respiration) is calculated using a combination of Farquhar kinetics [Farquhar et al., 1980] as modified by Collatz et al. [1991Collatz et al. [ , 1992, and the Ball-Berry equation (4) for stomatal conductance (g s , where g s = 1/r s ) [Ball, 1988] [see Sellers et al., 1996c, Appendix C].
where h s is relative humidity at the leaf surface and is related to relative humidity of the canopy air through r b and the transpiration rate, which is a function of wind speed and leaf geometry. Internal CO 2 concentrations (C s , C i , and C c ) are determined in an iterative loop that solves for internally consistent rates of A n and C c .
Mesophyll conductance (g m , where g m = 1/r m ) is calculated using the following equation: where vmax0 is the maximum rate of photosynthesis of sunlit leaves at the top of the canopy and is derived from the literature. Å is an integrating factor that expresses photosynthetic rate of the entire canopy as a function of the top leaf rate and the ability of the canopy to use photosynthetically active radiation (PAR). Å is strongly influenced by changes in NDVI. RSTFAC is a soil water stress factor  (see also Evans and Loretto [2000] and references therein for a discussion of the effects of soil water stress on mesophyll conductance). Four thousand is a constant used to adjust mesophyll conductance to achieve a drop in CO 2 pressure between C i and C c of about 8 Pa at high rates of photosynthesis [Evans and Loretto, 2000].
[9] C4 photosynthesis also discriminates against 13 C, but to a much lesser extent. In C4 plants, CO 2 is ''captured'' in the stomatal cavity and transported to the site of enzymatic fixation with PEPC. Since nearly all of the CO 2 that reaches the site of fixation in C4 plants is assimilated, the kinetic isotope effects that are expressed are largely those involved in transport. Consequently, in this model we assume that carbon isotopic discrimination in C4 plants is constant and equal to the isotope effect associated with diffusion through the stomatal pore, i.e., DPS C4 = 4.4%.
[10] There have been several previous studies modeling global terrestrial discrimination [Lloyd and Farquhar, 1994;Fung et al., 1997;Kaplan et al., 2002]. These studies differ from ours in the following ways (Table 1). Lloyd and Farquhar calculate monthly net assimilation as a simple function of climatologic monthly mean 0.5°Â 0.5°temperature and precipitation [Friedlingstein et al., 1992;Leemans and Cramer, 1991] and monthly 2.5°Â 2.5°wet and dry bulb temperatures from ECMWF. Landcover is from Wilson and Henderson-Sellers [1985]. Stomatal CO 2 concentrations are calculated using an analytic model of optimal stomatal behavior with respect to water use efficiency. Chloroplast CO 2 concentrations are based on observed relationships between canopy CO 2 , stomatal CO 2 , and chloroplast CO 2 . C3 discrimination is calculated using an equation similar to (3), although they include effects of autotrophic respiration and photorespiration and neglect the effect of transport across the leaf boundary layer. The effect of both of these differences on net discrimination is small. C3/C4 distributions in the work of Lloyd and Farquhar are based on a combination of physiological modeling and landcover maps from Wilson and Henderson-Sellers [1985]. C4 plants account for 21% of global GPP. C4 discrimination is a function of internal CO 2 concentrations and the leakiness of the bundle sheath cells. Net global C4 discrimination in the work of Lloyd and Farquhar is 3.6%. Total discrimination (C3 and C4) in the work of Lloyd and Farquhar is 14.8%; C3 discrimination is 17.8%.  Fung et al. [1997] calculate net assimilation at a 4°Â 5°resolution using monthly mean output from an older version of SiB2 coupled to a global climate model Randall et al., 1996]. Fung et al. calculate net assimilation-weighted C i /C a ratios by simultaneously solving for C i , stomatal conductance, and net assimilation. Discrimination is based on a two-step model of photosynthesis that includes only transport of CO 2 into the leaf, followed by enzymatic fixation. Fung et al. assume that savannas are 75% C4 and land areas covered by shrubs and groundcover are 50% C4. C4 plants make up 27% of GPP. C4 discrimination is fixed at 4.4%. Total discrimination in the work of Fung et al. is 15.7%; C3 discrimination is 20.0%.
[12] Kaplan et al. [2002] calculate net assimilation using a model of the terrestrial biosphere (BIOME4) driven by 0.5°Â 0.5°gridded climatology (Climate 2.2), which is an updated version of work by Leemans and Cramer [1991]. Vegetation, including C3/C4 distributions, is determined in a dynamic vegetation model with an agricultural mask. Kaplan et al. calculate water-limited C i /C a ratios in a process-based model of canopy conductance. The discrimination model is similar to that of Lloyd and Farquhar. C4 plant distributions in BIOME4 are determined in a dynamic vegetation model with an agricultural mask. C4 plants constitute only about 15% of global GPP, significantly less than suggested by the models here (J. O. Kaplan, personal communication, 2004), which range from 21% to 27% or the 18% estimated by Ehleringer et al. [1997]. C4 discrimination in the work of Kaplan et al. is determined in a manner similar to that of Lloyd and Farquhar. Discrimination for C4 plants is not given. Total discrimination in the work of Kaplan et al. is 18.6% or 18.1% with the agricultural mask; C3 discrimination is 20.0 ± 1.0%. The higher total discrimination in the work of Kaplan et al. relative to the other models, as well as observations [Bakwin et al., 1998;Miller et al., 2003], is probably due largely to differences in global C3/C4 distributions.
[13] Our simulation calculates net assimilation each 10min time step using an updated version of SiB2 coupled to assimilated 6-hourly 1°Â 1°weather from ECMWF for 1983 through 1993. We calculate internal CO 2 concentrations each time step in SiB2 as described earlier. One difference between our simulation and that of Fung et al. [1997] is that we calculate net assimilation solving for C c instead of C i . Discrimination is determined using equation (3). As Fung et al. point out, all other things being equal, neglecting isotope effects during transport of CO 2 from the stomatal cavity to the chloroplast can reduce total C3 discrimination by about 2 -3%. However, since we calculate net assimilation using C c instead of C i , net assimilation rates are lower and consequently the change in discrimination is much less than 2 -3%. We use C4 maps based on physiological modeling, satellite observations, and agricultural maps . C4 plants constitute 25% of GPP. C4 discrimination is fixed at 4.4%. Total discrimination in our study is 15.9%; C3 discrimination is 19.2%. A final feature of our discrimination module is that it can be coupled to a prognostic canopy air space and can thus be used to examine systematic variations between discrimination and transport on the time step of the simulation including recycling of respired CO 2 and Rayleigh-type distillation of canopy CO 2 during daytime photosynthesis.
[14] The d 13 C values reported here are weighted by assimilation. 12 C and 13 C fluxes are calculated separately for C3 and C4 plants and then summed to get grid values for each time step. At the end of the year, we divide by total net assimilation to retrieve the d 13 C of plant carbon assimilated within that grid cell during the preceding year. At the end of the 11-year run, we calculate mean and standard deviation of annual assimilation weighted d 13 C ratios for each 1°Â 1°g rid cell. Since we want to focus here on the response of the terrestrial discrimination to climate variations, we hold the concentration and carbon isotope ratio of canopy CO 2 constant. This necessarily neglects the influence of variations in d 13 C of canopy CO 2 on d 13 C of assimilated carbon due to (1) recycling of respired CO 2 , (2) enrichment of canopy CO 2 in 13 C due to preferential of removal of 12 C during photosynthesis, (3) turbulent mixing with the atmospheric boundary layer, and (4) seasonal variations in d 13 C of atmospheric CO 2 . A subsequent analysis will evaluate the influence of a dynamic canopy on d 13 C of terrestrial fluxes.

Spatial Variations in D 13 C of Plant Carbon
[15] Spatial variations in d 13 C of plant carbon are complex and reflect climate, water availability, biome type and C3/C4 distributions (Figure 4a). In general, the lightest d 13 C values for are found at high northern latitudes and reflect both the absence of C4 plants and increased discrimination in C3 plants due to zonal variations in the modeled C c /C a ratio. The starkest contrasts globally are those produced by variations in C3/C4 distributions. There is a zonal trend to this distribution and its influence on d 13 C; that is, the heaviest d 13 C values are found approximately 10°to 20°north and south of the equator, mirroring the distribution of C4 grasslands.  (1983)(1984)(1985)(1986)(1987)(1988)(1989)(1990)(1991)(1992)(1993) and (b) standard deviation in annual d 13 C values. During each time step, net assimilation is multiplied by the d 13 C of assimilated plant carbon. At the end of the year, the sum of net assimilation times d 13 C is divided by total annual net assimilation to give the d 13 C value of carbon assimilated in each grid cell during that calendar year. The mean is the average of the 11 annual values. The standard deviation is also for those 11 values. Please note the higher resolution used in reporting more negative mean values (Figure 4a). This is done in order to highlight spatial variations in C3 discrimination. There is tremendous spatial variability in d 13 C values of terrestrial carbon. The sharpest spatial contrasts are the result of differences in C3/C4 ratios. Variation in d 13 C values of C3 plants is largely driven by differences in relative humidity during the growing season. Larger than average standard deviations are found in areas where ECMWF predicts significant variations in annual precipitation. See color version of this figure at back of this issue.
However, there are also important longitudinal patterns in d 13 C values that are produced not only by C3/C4 distributions, but also by variations in climate and biome type. Examples of this include the influence of high C4 ratios in the North American Plains and African Sahel on d 13 C of terrestrial biomass of the Western and Eastern Hemispheres. Table 2 presents d 13 C ratios of all plants (d 13 C Tot ) and C3 plants only (d 13 C C3 ) broken down by biome.
[16] The standard deviations in annual d 13 C values for individual grid cells for the 11 years are generally less than 0.3% (Figure 4b). The largest standard deviations are generally found in areas where ECMWF predicts large variations in annual precipitation and occasional drought. The most prominent example is the western and southern edges of the Amazon basin, which in ECMWF reanalysis experienced significantly drier conditions in 1983 -1986, relative to the rest of the 11-year simulation.

C3 Discrimination and Vapor Pressure Deficit
[17] Laboratory and field studies show that discrimination in C3 plants is largely controlled by the influence of vapor pressure deficit (vpd) on stomatal conductance and C c /C a ratios [Farquhar and Richards, 1984;Farquhar et al., 1982Farquhar et al., , 1989Berry, 1988;Condon et al., 1993;Hall et al., 1993]. Since vpd of canopy air changes over the course of the day and from one day to the next, carbon isotope discrimination should change just as quickly. Until recently, however, it has been difficult to quantify changes in discrimination in natural ecosystems on these short timescales.
[18] Nighttime Keeling plot intercepts record d 13 C of respired CO 2 fluxes (d 13 C R ) from aboveground (bole, stem, leaf, and heterotrophic) and belowground (root and decomposing soil organic matter) sources [Keeling, 1958]. Högberg et al. [2001] have shown that tree girdling reduced soil respiration by up to 37% after only 5 days, indicating that roots are responsible for approximately 40% of the total soil CO 2 flux. Using the fact that d 13 C of root respiration reflects the carbon isotope ratio of recent photosynthate [Andrews et al., 1999;Lin et al., 1999], Ekblad and Högberg [2001] and Bowling et al. [2002] then demonstrated that variations in d 13 C of ecosystem respiration are strongly correlated with changes in relative humidity and/or vapor pressure deficit that have occurred within the previous 2 to 10 days. Consequently, changes in nighttime Keeling plot intercepts reflect humidity-induced changes in discrimination through the release of recently fixed carbon during autotrophic respiration.
[19] In order to evaluate the ability of our model to accurately represent the relationship between C3 discrimination and relative humidity, we compare simulated mean daily vapor pressure deficits and assimilationweighted d 13 C values for C3 plants to the observed relationships ( Figure 5). We have chosen sites similar to those presented by Ekblad and Högberg [2001] and Bowling et al. [2002], as well as others spanning a range of vapor pressure deficits and biomes. In general, the model successfully captures the relationship between vpd and d 13 C and is bracketed by the observations. However, there are differences between model results and observation that are worth noting. The d 13 C values in the simulation are similar to those observed by Ekblad and Högberg, but slightly less negative than those of Bowling et al. There are a couple of reasons that it is difficult to assess the importance of these absolute differences in d 13 C values. First, both sets of observations include undetermined offsets due to the fact that the d 13 C of decomposing soil organic matter (SOM) is unknown. In general, CO 2 from decomposing SOM is thought to be enriched in 13 C relative to recently fixed photosynthate as a consequence of a terrestrial ''Suess Effect'' [e.g., Ciais et al., 1999]. On the basis of flux-weighted ages of decomposing SOM [Fung et al., 1997] and a long-term record of d 13 C of atmospheric CO 2 , we might expect d 13 C of SOM in northern Sweden to be enriched relative to recent photosynthate by about 0.7%, whereas the sites of Bowling et al. would be slightly less enriched and more variable (0.2% to 0.7%). However, these are only rough estimates, and more would have to be known about the age and turnover rate of the SOM at specific sites to be more precise. If decomposing soil organic matter produced 60% of the soil CO 2 flux, a 0.7% offset would mean that d 13 C of recently fixed carbon was approximately 0.4% more negative than observations indicate. Another problem with comparing the two studies is that they use slightly different methodologies. The differences include: (1) different methods were used for collecting and analyzing gas samples, (2) Ekblad and Högberg looked at a single site, whereas Bowling et al. used   deficit in the canopy air is a dominant control on discrimination on diurnal timescales.

C3 Discrimination and Monthly Precipitation
[20] Ometto et al. [2002] examine the relationship between C3 discrimination and monthly precipitation in tropical forest sites of the Amazon basin. In data collected near Santarém, Brazil, there is an inverse relationship between d 13 C of ecosystem respiration (d 13 C R ) and monthly precipitation that persists at rain rates up to 300 mm/month. At much higher rates of precipitation, i.e., greater than 400 mm/month, the relationship appears to reverse. They argue that at the lower rain rates, changes in d 13 C R are controlled by stomatal-induced variations in C3 discrimination caused by a combination of soil water stress and relative humidity. It is unclear why discrimination decreases at precipitation rates greater than 400 mm/month. Ometto et al. also collected data from a second tropical forest site near Manaus, Brazil. At Manaus, the relationship between d 13 C R and monthly precipitation is less clear, even though the site is similar in many respects to that of Santarém.
[21] We examined the relationship between monthly precipitation and C3 discrimination in the simulation for the grid cell containing Santarém, Brazil, as well as a temperate forest grid cell containing the WLEF tall tower site in Wisconsin ( Figure 6). In the simulation, C3 discrimination increases with increasing monthly precipitation. The relationship is slightly nonlinear, with a steeper slope at very low rates of precipitation. In other words, changes in Figure 5. Relationship between vapor pressure deficit (vpd), i.e., the difference between saturated water vapor pressure and vapor pressure of the canopy air, and d 13 C of ecosystem respiration (d 13 C R ) by C3 plants. The two lines are from observations. The data points are from the simulation. The solid line is the proposed linear fit to measurements from a mixed coniferous boreal forest in northern Sweden (y = 0.2242x À 27.3; r2 = 0.94, n = 5) [see Ekblad and Högberg, 2001, Figure 3]. Meteorologic data were collected at a weather station 200 m from the site. The vpds are calculated from the 2-day mean for maximum air temperature and the 2-day mean relative humidity at the maximum air temperature. The dotted line is a simple logarithmic fit to data proposed by Bowling et al. [2002] (y = 2.54Lnx À 27.8; r 2 = 0.79, n = 22; see Figure 7). The curve is constructed from 22 points from four sites in Oregon. Bowling et al. calculated vpd using either on-site instrumentation or data from Oregon Climate service stations as near to the site as possible. The d 13 C R values are constructed from 2-day averages of Keeling plot intercepts using geometric mean regressions of nocturnal data. Since there is likely to be a delay between the time the carbon is photosynthetically fixed and the time that it is respired, in the observations a time lag was selected by looking for the best statistical fit between d 13 C R and 2-day mean vpd over varying time periods. Although Ekblad and Högberg used a 3-to 4-day time lag for the relationship between d 13 C R and relative humidity, a 2-to 3-day lag was best for fitting d 13 C R and vpd; Bowling et al. used time lags ranging from 5 to 10 days, depending on the site. Data points are from the simulation and are mean daily vpd for daylight hours only and assimilation-weighted d 13 C values. The red solid circles, blue stars, and green crosses are from grid boxes in Oregon which are similar to the locations in the Bowling et al. study. Other data points are from sites and months intended to represent a range of conditions from very wet to very dry. The locations are as follows: WLEF (45.5°N, 90.5°W); W Wyoming (41.5°N, 109.5°W); Manaus (2.5°S, 55.5°W); N Oregon (48.5°N, 124.5°W); and S Oregon (43.5°N, 123.5°W). In the simulation, d 13 C of canopy CO 2 is held constant at À7.9%.
precipitation rates have a greater impact on carbon isotope discrimination under very dry conditions than they do when the environment is already wet. Discrimination at WLEF is about 1% less than that predicted for Santarém, largely because of the greater relative humidity of the tropical forested site. Model results agree with the observations to the extent that we predict a positive correlation between monthly precipitation and discrimination; however, there are several important differences between the model and observations. First, observed d 13 C R for both Santarém and Manaus is depleted in 13 C by up to 2% relative to modeled discrimination. Second, the slope of the relationship in the observations at Santarém is about 3 times greater than the slope predicted for the grid cell containing Santarém. In the observations, d 13 C of ecosystem respiration decreases by approximately 3% for changes in precipitation from 50 to 300 mm/month; in the simulation, the change is closer to 1%. Third, at Santarém, d 13 C of ecosystem respiration returns to heavier values at precipitation rates greater than 400 mm/month, whereas modeled discrimination continues to increase with increasing precipitation rates.
[22] Reasons for the enrichment in modeled values relative to observations at both Santarém and Manaus will be addressed later in the paper. With respect to the relationship between discrimination and monthly precipitation, we think that there may be several reasons for the differences between the model and observations at Santarém. First, while we have no reason to doubt the results of Ometto et al. [2002], the question must be asked: Is this site representative of tropical forest behavior in general? Data collected at Manaus do not refute the relationship observed at Santarém; however, neither do they confirm it. In part, this is because of the narrow range of monthly precipitation rates observed at Manaus; however, more data will have to be collected before we can say that the behavior at Santarém is characteristic of tropical forests in general. A second possibility is that there are important natural variations accompanying changes in monthly precipitation that are not captured in the model. The two most probable candidates are relative humidity and/or soil water. Since SiB2 seems to do a good job simulating the relationship between discrimination and vapor pressure deficit, the problem may be with the way that we model soil water stress. Throughout the 11 years of the simulation, there is no significant soil water stress on plant growth in the grid cell containing Santarém, even during the driest months. Unfortunately, we have no observations to directly compare to the simulation. Terrestrial ecosystem models sometimes predict decreases in photosynthetic rates during the tropical forest dry season due to drought stress; however, this decline in assimilation rates is not always supported by the observations [Tian et al., 1998;Saleska et al., 2003]. Finally, there is a similar Figure 6. Relationship between monthly precipitation and C3 discrimination. Discrimination has been converted from del values using the formula: Discrimination = (d 13 C air À d 13 C R )/(1 À d 13 C air /1000); where d 13 C air = À7.9% and d 13 C R is the carbon isotope ratio of ecosystem respiration. Triangles and asterisks are from measurements in primary evergreen tropical forests: triangles (Santarém, Brazil), asterisks (Manaus, Brazil) [see Ometto et al., 2002, Figure 7]. The error bars represent the standard error of the isotope measurement. The curve is a second-order polynomial fit to the data from Santarém. The open circles are simulated discrimination for the grid box containing Santarém. Simulated data from Manaus are nearly identical. The solid squares are simulated discrimination for summer months for the grid box containing WLEF tower in Wisconsin. The lines through the simulated data are logarithmic fits to the data.  [Pataki et al., 2003]. Error bars are standard deviations. SiB2's broadleaf evergreen biome (À26.6 ± 0.3%, n = 1098) is compared to observations from tropical forests (À27.6 ± 0.94%, n = 5); broadleaf deciduous (À25.9 ± 1.2%, n = 318) is compared to temperate broadleaf (À25.1 ± 1.8%, n = 18); broadleaf and needleleaf (À26.3 ± 0.4%, n = 769) is compared to temperate conifer (À26.6 ± 1.2%, n = 13); and needleleaf (À27.1 ± 0.4%, n = 2017) and needleleaf deciduous (À26.6 ± 0.4%, n = 952) are combined (26.9 ± 0.5%, n = 2969) and compared to boreal (À26.3 ± 0.6%, n = 7). In general, simulated d 13 C values are close to observed values although there is slightly less variability in simulated values than in observations. In addition, simulated d 13 C values in the tropics are heavy relative to observed values. pattern between observed mean carbon isotope composition of ecosystem respiration and mean annual precipitation; that is, discrimination increases with increasing precipitation up to a point, where it then returns to lower values at very high rain rates. In the simulation, annually integrated d 13 C values are generally more negative in wet climates, and they do not return to less negative values in extremely wet climates.

C3 Discrimination by Biome
[23] The number of measurements of d 13 C R that have been collected is constantly increasing and can serve as an additional method for evaluating simulation results. Data for Figure 7 are taken from Pataki et al. [2003] and show mean d 13 C R values and standard deviations for four different biomes: Broadleaf, Temperate Broadleaf, Temperate Conifer, and Boreal. We have plotted mean d 13 C R values and standard deviations for comparable biomes from our simulation against the observed values. In general, simulated mean values are relatively close to the observed means, and the standard deviations encompass the range of observed values. However, there are a couple of important differences between the simulation and the observations. First, the mean model d 13 C R for the tropical forests is 1% less negative than the mean of the observations. If simulated discrimination in the tropical forests is incorrect, it could be the result of a couple of factors. First, our estimate for maximum C3 discrimination (28.2%) may be too small. This could be because the isotope effect that we use for CO 2 fixation with Rubisco (30%) is not large enough, or because of our assumption that 7.5% of the carbon fixed by C3 plants occurs through reactions with PEPC is too great. Either one is a possibility. Measurements of isotope effects associated with rubisco are precise, but the values range from approximately 26% to 38% [O'Leary, 1981;Brugnoli and Farquhar, 2000]. With respect to the amount of CO 2 in C3 plants that reacts with PEPC, very little is known. A maximum of 10% has been suggested [Raven and Farquhar, 1990] and Lloyd and Farquhar [1994] assumed that it was 5%, so perhaps our estimate of 7.5% is too high. Changing from 7.5% to 5% reduces discrimination by about 0.6%. Of course, while the change improves the match to observed d 13 C R for tropical forests, it hurts the fit with Boreal and Temperate Broadleaf biomes. Another possibility is that there are biome-specific differences in this fraction of which we are unaware. A second reason why our d 13 C values for tropical forests are enriched in 13 C relative to observations could be due to the fact that we do not include the effects of recycling of isotopically depleted respired CO 2 . In this study we hold d 13 C of canopy CO 2 constant at À7.9%. Estimates of the amount of respired CO 2 that is recycled are as high as 30 to 40% [Sternberg et al., 1997], although others have argued that it is actually quite negligible [Lloyd et al., 1996[Lloyd et al., , 1997. In sensitivity tests in which we coupled the biosphere to a dynamic canopy, recycling was greatest in tropical forests, which is consistent with observations, and caused d 13 C of assimilated carbon to be depleted in 13 C by 1-4%. Working against this effect, however, is an offset due to the terrestrial ''Suess Effect.'' On the basis of flux-weighted ages of decomposing SOM [Fung et al., 1997] and a long-term record of d 13 C of atmospheric CO 2 , the terrestrial Suess Effect in tropics is estimated to be about +0.6%, which would counter some of the effects of recycling. A second and perhaps more important difference between the observations and the model is that there appears to be a wider range of values in the observations than in the model. Unfortunately, there are still so few observations that it is difficult to know how representative they are. Nonetheless, it is entirely possible that the wider range of values in the observations is due to heterogeneity within biomes that is not captured by the model.

Zonal Trends in D 13 C of Plant Carbon Compared to Atmospheric Data
[24] Results of the simulation can also be compared to atmospheric data. Miller et al. [2003] have analyzed highfrequency deviations from a seasonal spline fit to time series of CO 2 and d 13 C measurements from continental stations from the NOAA Globalview Flask Network in order to deduce the carbon isotopic ratios of growing season CO 2 fluxes in the Northern Hemisphere. They argue that these deviations from the seasonal signal are produced by fluxes to and from the local vegetation and are therefore a record of the isotopic ratio of that vegetation. Miller et al. have plotted the d 13 C values from these sites as a function of latitude, as well as zonal mean d 13 C values from Fung et al. [1997] and Lloyd and Farquhar [1994]. Their figure is recreated here with the addition of our zonal mean assimilation-weighted d 13 C values for the Northern Hemisphere growing season (Figure 8). While it is difficult to make a direct comparison, because the observations largely represent regional values of terrestrial discrimination, it appears that our zonal mean is an improvement over previous simulations in that it captures some of the flatness of the signal between 20°N and 40°N, as well as the rapid descent to more negative values north of 40°N. Furthermore, although we are still limited by a paucity of data, both observations and simulated results indicate that spatial distribution of d 13 C values of the terrestrial biosphere are more complex than previously thought.

Conclusion
[25] A multistage model of carbon isotope discrimination during photosynthesis and global maps of C3/C4 plant ratios coupled to an ecophysiological model of the terrestrial biosphere driven by observed meteorology and constrained by satellite-derived NDVI has been used to predict the carbon isotope ratios of terrestrial plant carbon. Simulated mean annual C3 discrimination for 1983-1993 is 19.2%; total mean annual discrimination by the terrestrial biosphere (C3 and C4 plants) is 15.9%.
[26] It is very difficult to constrain the details of the simulated results because we are still sorely lacking in observations. However, within these limitations we find the following. The modeled response to changes in vapor pressure deficit compares favorably to observations. Simulated monthly discrimination in tropical forests is less sensitive to changes in precipitation than is suggested by some observation. Mean model d 13 C ratios for Broadleaf, Temperate Broadleaf, Temperate Conifer, and Boreal biomes compare well with the few measurements available; however, there is more variability in observations than in the simulation, and modeled d 13 C values for tropical forests appear heavy relative to observations. Finally, simulated zonal d 13 C ratios in the Northern Hemisphere capture the complexity of the zonal d 13 C inferred from atmospheric measurements better than previous investigations.
[27] Acknowledgments. We would like to thank Inez Fung for helpful discussions. This work was funded by the NASA Earth Science Enterprise Interdisciplinary Science Program under contract NAS-31730 and the National Science Foundation under contract 00223464. Figure 2. Distribution of C3 and C4 plants  is determined using remote sensing products [Defries and Townshend, 1994;DeFries et al., 1998DeFries et al., , 1999aDeFries et al., , 1999bDeFries et al., , 2000, physiological modeling [Collatz et al., 1998], crop fraction maps [Ramankutty and Foley, 1998], and areal coverage data on crop types [FAO, 1998].  [1983][1984][1985][1986][1987][1988][1989][1990][1991][1992][1993] and (b) standard deviation in annual d 13 C values. During each time step, net assimilation is multiplied by the d 13 C of assimilated plant carbon. At the end of the year, the sum of net assimilation times d 13 C is divided by total annual net assimilation to give the d 13 C value of carbon assimilated in each grid cell during that calendar year. The mean is the average of the 11 annual values. The standard deviation is also for those 11 values. Please note the higher resolution used in reporting more negative mean values (Figure 4a). This is done in order to highlight spatial variations in C3 discrimination. There is tremendous spatial variability in d 13 C values of terrestrial carbon. The sharpest spatial contrasts are the result of differences in C3/C4 ratios. Variation in d 13 C values of C3 plants is largely driven by differences in relative humidity during the growing season. Larger than average standard deviations are found in areas where ECMWF predicts significant variations in annual precipitation.