Role of methane and biogenic volatile organic compound sources in late glacial and Holocene fluctuations of atmospheric methane concentrations

2.1. Choice of Time Slices

[11] Although the surface fluxes produced by BIOME4-TG were available at every 1 ka time slice from 21 ka to 0 ka, computational requirements limited our selection of time slices for the LMDz-INCA simulations to four key periods in the past: LGM (21 ka), early Holocene (10 ka), mid-Holocene (6 ka), and preindustrial conditions (0 ka). This preindustrial control simulation represents potential natural conditions, without the impact of anthropogenic land use or emissions from fossil fuel use. The selected paleoclimate time slices for the chemistry-climate model experiments were chosen because these periods are the focus of modeling efforts in the Paleoclimate Model Intercomparison Project II (PMIP II) [Harrison et al., 2002], and because they represent extremes in climate forcing.

[12] At the LGM, the Earth's orbital configuration was similar to today but the continental ice sheets were at a maximum and greenhouse gas concentrations were at a minimum [Dyke and Prest, 1987; Raynaud et al., 1993; Svendsen et al., 2004]. In addition to the large changes in terrestrial geography, the ocean surface was colder and the distribution of sea ice was greatly expanded compared to the present [CLIMAP Project Members, 1976, 1981]. The configuration of the Earth's orbit was substantially different from today or the LGM during the intervening period of the early to mid-Holocene. The phasing of the precession (19 kyr, 23 kyr) and obliquity (41 kyr) cycles was such that at the early Holocene (circa 11 ka) the high latitudes of the Northern Hemisphere received more insolation than any other time over the past 21 ka, both during summer and annually [Berger et al., 1993]. At 10 ka, the northern high-latitude insolation was near its maximum; at 6 ka this effect was still significant but somewhat smaller. As a direct consequence of the orbital forcing, parts of the northern high latitudes experienced summers in the early and mid-Holocene that were considerably warmer than the present [Bigelow et al., 2003; MacDonald et al., 2000; Ritchie et al., 1983].

[13] Other boundary conditions were different between 10 ka and 6 ka: the Laurentide ice sheet, although substantially reduced from its LGM size, was still sufficiently large to have a major downwind cooling effect during the early Holocene [Harrison et al., 1992; Mitchell et al., 1988]. At 10 ka, a residual Fennoscandian ice sheet also persisted in northern Europe, and sea level was ∼30 m lower compared to 6 ka. Thus northern Europe and eastern North America experienced a thermal maximum several thousand years after the insolation maximum. In contrast, atmospheric greenhouse gases (CO₂, CH₄, and N₂O) had reached near preindustrial concentrations by 10 ka [Raynaud et al., 1993], though they continued to fluctuate through the Holocene [Indermühle et al., 1999].

2.2. BIOME4-TG Vegetation Model

[14] BIOME4-TG is an equilibrium-state terrestrial biosphere model that couples biogeography and biogeochemistry and is able to calculate the vegetation distribution at any time in the past or future in response to both changing climate and atmospheric CO₂ concentrations. The model is driven by monthly climatological mean temperature, precipitation, and solar radiation [Kaplan, 2001]. BIOME4-TG simulates both vegetation structure in the form of 28 world biomes, and a suite of standard biogeochemical variables, including NPP, LAI, autotrophic and heterotrophic respiration, water runoff, and soil moisture. The biomes in BIOME4-TG are defined so as to be downward compatible with the classification used in the BIOME 6000 data set of reconstructed past vegetation [Prentice et al., 2000; Prentice and Webb, 1998]. BIOME4 has been used extensively for simulating past vegetation distributions and evaluating paleoclimate simulations by GCMs [Harrison, 2000; Kaplan et al., 2003].

[15] BIOME4-TG extends the functionality of BIOME4 by including modules to simulate wetland area, CH₄ emissions from wetlands, and BVOC emissions from all vegetation types. To determine wetland area, BIOME4-TG uses simple algorithms based on topography and soil moisture; CH₄ emissions are based on ecosystem carbon turnover in wet soils [Christensen et al., 1996; Kaplan, 2002]. Biogenic VOC emissions are calculated in BIOME4-TG on the basis of algorithms from the Global Model of Natural Volatile Organic Compound Emissions [Guenther et al., 1995].

2.3. Climate Scenarios Used by BIOME4-TG

[16] Paleoclimate scenarios used by BIOME4-TG were prepared by combining a GCM paleoclimate simulation with a baseline climatology of 20th century mean climate. Paleoclimate simulations were performed with the Hadley Centre Unified Model (HC-UM), including a model of the ocean mixed layer (“slab ocean”) [Kaplan et al., 2002]. The paleoclimate scenarios were produced by applying climate anomalies, that is, differences in monthly mean climates between the HC-UM paleoclimate and control (“present-day”) simulations, to a “baseline” representing late 20th century mean climate. The HC-UM simulations were performed at 1000-year intervals from 21 ka to “0 ka” (preindustrial time). Boundary conditions for the climate simulations were orbital forcing, physiography (orography, ice sheets, sea level), and atmospheric CO₂ concentration. Land surface conditions (vegetation type, density, height, roughness, and phenology), including potential vegetation feedbacks to the climate system, were calculated through asynchronous coupling to the BIOME4 model [see, e.g., de Noblet et al., 1996; Kaplan, 2006]. Because of the 1000-year time slice approach of the GCM simulations, rapid climate change events, for example, the Younger Dryas and 8.2 ka event, could not be taken into account in the paleoclimate scenarios. While these transient phenomena are of great interest in interpreting the ice core record of atmospheric CH₄, we limit our analysis to the longer-term changes in the CH₄ cycle, which we hypothesize to be affected by different processes.

[17] BIOME4-TG was run on a 0.5° × 0.5° geographic grid, corresponding to the resolution of the baseline climatology. The CLIMATE 2.2 baseline climatology is a gridded long-term mean climatology for the 20th century that includes monthly mean temperatures and percent of potential sunshine hours, and monthly total precipitation [Kaplan et al., 2003; Leemans and Cramer, 1991]. The grid includes all “virtual” land grid cells on the continental shelf areas that were exposed at the LGM and times between then and the present. The land and ice mask at each 1 ka interval was defined by a series of simulations from the dynamic ice sheet–crustal deformation model ICE-3G [Peltier, 1994]. The ICE-3G output was linearly interpolated from its native 1° grid to the 0.5° grid used by BIOME4-TG.

[18] An atmospheric CO₂ concentration of 324 ppm was used to force BIOME4-TG for the present-day simulation. This value represents the mean [CO₂] during the period of measurement of the station data upon which the climatology is based. For paleoclimate simulations, CO₂ concentrations were interpolated at 1 ka intervals from the Taylor Dome and Byrd Antarctic ice core CO₂ records [Neftel et al., 1988; Indermühle et al., 1999; Stauffer et al., 1998].

2.4. LMDz-INCA Atmospheric Chemistry-Climate Model

[19] LMDz (Laboratoire de Météorologie Dynamique, zoom) is a full complexity, 3D General Circulation Model (GCM) initially developed for climate studies by Sadourny and Laval [1984]. As part of the IPSL Earth System Model, the GCM has been applied in land surface–climate feedback studies by Friedlingstein et al. [2001] and Dufresne et al. [2002]. In LMDz the finite volume transport scheme of Van Leer [1977] is used to calculate large-scale advection of tracers and deep convection is parameterized according to the scheme of Tiedke [1989]; a local second-order closure formalism has been integrated to describe turbulent mixing in the planetary boundary layer (PBL). LMDz (version 3.3) has a horizontal resolution of 3.8 × 2.5 degrees on the horizontal and uses 19 vertical σp levels extending from the surface to 3 hPa.

[20] The Interactive Chemistry and Aerosols (INCA) model has been integrated into LMDz. INCA includes 83 chemical species, 43 photolytic reactions, and 217 thermochemical reactions. INCA has been designed to simulate tropospheric chemistry, emissions, and deposition of primary trace species of the troposphere including non-CH₄ hydrocarbons (NMHC) and non-CH₄ volatile organic compounds (NMVOC). Detailed description and evaluation of LMDz-INCA and emission inventories can be found in Hauglustaine et al. [2004], Folberth et al. [2005a, 2005b], and Lathière et al. [2005].

2.5. Emission Inventories

[21] In order to represent preindustrial conditions, the surface emission magnitude for several species in the standard LMDz-INCA emission inventory, which is characteristic for the present-day environment, had to be readjusted. Modifications follow the representative inventory for the year 1850 of Lelieveld et al. [1998] for most species, with Yienger and Levy [1995] for NO_x. The changes include the elimination of all surface emissions from fossil fuel combustion, a reduction of the biomass burning emissions (including agricultural waste burning and fuel wood use) for all species to 15% of their present-day magnitude, and a decrease in the NO_x soil emissions to 2.9 Tg N yr⁻¹ (as opposed to 5.5 Tg N yr⁻¹ for the present). The global atmospheric NO_x source from lightning was kept at its present-day value of 5.0 Tg N yr⁻¹. CH₄ emissions from rice paddies, wild animals and ruminants, and landfills (not included in BIOME4-TG) were set to 33% (19.7 Tg C yr⁻¹), 29% (20.3 Tg C yr⁻¹), and 25% (6.7 Tg C yr⁻¹), respectively, of their present-day magnitude, following Lelieveld et al. [1998].

[22] BIOME4-TG simulated inventories of wetland CH₄ and biogenic isoprene, terpenes, and other reactive volatile organic compounds (ORVOC). Because LMDz-INCA requires surface fluxes of individual VOC species as input, the geographic distribution of all species other than isoprene and terpenes was determined by the simulated distribution of total ORVOC. Emission magnitude was determined by rescaling the standard emission inventory for each species by the total ORVOC flux simulated by BIOME4-TG at each time slice. In other words, we assume a linear relationship between total ORVOC emissions and the emissions of any given species; an assumption that, in view of the lack of any other information, seems reasonable. Table 1 summarizes the emission inventories for the key species as used in this study.

[23] In the three paleo time slice experiments it was assumed that the remaining anthropogenic portion of preindustrial surface RVOC emissions did not yet exist. This includes CH₄ emissions from rice paddies and landfills, and CO emissions from agricultural waste burning and fuel wood use. CH₄ emissions from ruminants were kept at their preindustrial magnitude. The elimination of anthropogenic emissions amounts to approximately 15% and 30% of the total preindustrial fluxes of CH₄ and CO respectively, but in each case constitute only slightly more than 2% of the total preindustrial RCC emission magnitude of ∼1110 TgC yr⁻¹. We do not consider our emission inventories to be an exact representation of the conditions at the four time periods considered here, but rather a rough estimate based on many assumptions, which nevertheless seem reasonable in view of the present knowledge and the purpose of this study.

2.6. Experiments Performed by LMDz-INCA

[24] At the four time slices described above, we ran LMDz-INCA to calculate OH burdens, and CH₄ lifetime and concentrations. While the emission inventories of CH₄ and BVOC simulated by BIOME4-TG varied geographically with changing sea level and ice sheet conditions, all other forcings on LMDz-INCA were held constant in order to specifically assess the role played by changing surface emissions on atmospheric CH₄ concentrations.

[25] The boundary conditions affecting the physical climate, including insolation and land and sea ice distributions, were fixed to preindustrial conditions in each LMDz-INCA experiment. This resulted in similar temperatures, winds, water vapor, and surface albedo simulated by the chemistry-climate model. The stratospheric ozone column was also kept constant between the time slices.

[26] Each simulation was initialized by a set of restart conditions from a previous model run and was taken over a 30 months period to allow a new equilibrium state in reactive trace gas concentrations to be reached. To facilitate an estimation of the photochemical equilibrium for CH₄ – this would otherwise take approximately 20 model years – the CH₄ field was globally initialized to a value close to the expected equilibrium concentration. The global mean CH₄ concentration was closely monitored and found to be in equilibrium within reasonable limits (continuous increase or decrease in CH₄ of less than 5% from year 1 to year 2) in all four cases. The change in the global annual mean CH₄ concentration from year 1 to year 2 amounts to 0.5%, −0.1%, 1.1%, and 3.1% for the PI, 6 ka, 10 ka, and LGM run, respectively.

[27] Interannual variability in these runs is expected to be on the order of 0.5 to 1% for both global mean CH₄ concentration and OH burden. Hence a substantial amount of the change between year 1 and year 2 can be explained by the model's interannual climate variability, which was assessed by analyzing a four model–year record of control runs for present-day conditions. All other species have photochemical lifetimes sufficiently short to reach equilibrium well within six model months. The last 12 months of each run were selected for analysis.

[28] The boundary conditions and model setup described above have significant implications for the modeled global distribution of OH and consequently for CH₄ concentrations. Uncertainties resulting from this modeling approach are discussed fully in the results section.

3.1. Vegetation Distribution

[29] Qualitative comparison of global vegetation distribution simulated in our control “preindustrial” simulation (Figure 2a) with modern observations, descriptive accounts of presettlement vegetation, paleodata synthesis, and recent modeling studies yields good agreement in most parts of the world [Haxeltine and Prentice, 1996; Prentice and Webb, 1998; Walter, 1973]. Boreal, temperate, and tropical forests occupy areas near their observed limits along with deserts and tropical savannas. BIOME4-TG, similar to other global vegetation models, systematically fails to reproduce observed vegetation distribution at locations where it underestimates the area of grassland and temperate savannah in favor of forest, for example, in central North America and the South American pampas [Sitch et al., 2003]; the model also overestimates the area of forest in hypermaritime boreal regions [Kaplan et al., 2003]. These biomes may be the result of processes including long-term regular disturbance or extreme diurnal variability in meteorological conditions [Transeau, 1935; Walter et al., 1975], processes that are not incorporated in BIOME4-TG. On the long timescales investigated in this study, however, shifts in global vegetation were much larger than the uncertainty associated with the simulation of present-day distributions [Prentice, 2001]. Over the past 21 ka global vegetation distribution responded to the melting ice sheets, rising sea level, changing climate and increasing atmospheric CO₂ concentrations that accompanied deglaciation (Figures 2b, 2c, and 2d).

[30] Large ice sheets, tundra, and xeric herbaceous and scrub vegetation dominated the LGM landscape (Figure 2d). The LGM simulation shows the prevalence of graminoid and forb tundra in Siberia and on the Bering land bridge at the LGM that are observed in paleoecological reconstructions of vegetation during this time [Bigelow et al., 2003]. With no widespread modern analog, the graminoid and forb tundra is indicative of very cold and dry climate conditions that favor herbaceous vegetation over shrubs. As both paleodata synthesis and other modeling studies have suggested, significant areas of tropical rain forest persisted at the LGM, particularly in the most humid areas of the Sunda shelf, and Amazon and Congo basins [Cowling et al., 2001; Hope et al., 2004; Prentice et al., 2000]. Otherwise, the tropics are dominated by xerophytic scrub vegetation, which may have occurred as a consequence of cool temperatures, low CO₂, or a combination of the two [Jolly and Haxeltine, 1997].

[31] During the early and mid-Holocene (10 ka and 6 ka time slices), boreal regions are characterized by a large increase in the area of forest compared to LGM (Figures 2b and 2c). The Eurasian boreal forest reaches its farthest northern extent at 10 ka, coincident with the summertime insolation maximum at this time. Further south in central Siberia the same time period is marked by the persistence of graminoid and forb tundra in the model simulation. While perhaps not a classically defined tundra under these 10 ka conditions, this biome type is nonetheless indicative of the cold, xeric, and low-productivity environment. In North America, the persistence of the Laurentide ice sheet prevented the boreal forest from reaching its maximum northward extent until circa 6 ka. Tropical vegetation developed to conditions similar to the present already by the early Holocene, and at both 10 ka and 6 ka forests indicative of a more humid climate compared to preindustrial developed in South Asia because of the intensified summer monsoon. In contrast, the Sahara desert region does not show the shift in vegetation recorded in the paleorecord [Hoelzmann et al., 1998]; this may be because the dynamic interplay between Holocene sea surface temperatures, land surface conditions, and the recycling of precipitation necessary to produce the vegetation shift could not be fully captured by the HC-UM GCM simulations [Harrison et al., 1998].

3.2. Wetland Area and CH₄ Emissions

[32] The evolution of wetland areas from the LGM to the present is related to both the development of boreal peatlands on formerly glaciated areas, and the submergence of wetlands on both tropical and boreal continental shelves that were exposed at the LGM. While total global wetland area decreased slightly between the LGM and preindustrial, CH₄ emissions remained relatively constant over this time, with only a small maximum during the late Pleistocene (Figure 3). The change in global CH₄ emissions is caused by the simultaneous decrease in area of dominant tropical sources, and the overall increase in CH₄ emissions due to warmer and wetter global climate and higher CO₂ concentrations that stimulated vegetation productivity.

[33] In our simulations for the preindustrial, global wetlands mainly appear in tropical and boreal lowland regions, with their largest extents in the Amazon basin and western Siberia, respectively. Other major wetland areas occur in the Pantanal of southern South America, in boreal North America, Eastern Europe, and in the Congo basin and Indonesian lowlands (Figure 4a). Total global wetland area at the preindustrial is calculated to amount to 5.8 × 10⁶ km², which is within the range of other estimates, although these figures vary widely even for the present day [Kaplan, 2002]. BIOME4-TG probably underestimates the concentration of wetland area in boreal Canada and specifically in the Hudson Bay Lowlands (HBL). However, recent studies suggest that the HBL is a relatively small CH₄ source (∼0.5 Tg yr⁻¹) when integrated over the whole region [Roulet et al., 1994; Worthy et al., 2000]. We therefore expect that the underestimate of Canadian wetland area by BIOME4-TG does not affect global CH₄ emissions or their geographic distribution significantly.

[34] At the LGM, wetlands were greatly reduced from their current extent in North America, Europe and western Siberia because of the presence of the ice sheets and perennially frozen soils (Figure 4d). Large wetland areas are simulated on the exposed continental shelves, however, particularly in Beringia, and on the tropical Sunda and Gulf of Carpenteria shelves. Total global wetland area at the LGM is simulated at 6.8 × 10⁶ km², slightly more than at preindustrial, but within the uncertainty of present-day estimates.

[35] From the LGM to the early Holocene (circa 12 ka), wetland area remained relatively constant and then decreased slightly through the course of the Holocene. At 10 ka wetlands are seen occupying much of western Siberia and Europe near preindustrial extent, while some wetland area still remained on the exposed remnants of Beringia and the Sunda shelf, and most of the boreal eastern North America was under glacial ice (Figure 4c). At 6 ka, with the final melting of the Laurentide ice sheet and the stabilization of sea level, wetland distribution is very similar to the preindustrial (Figure 4b).

[36] In contrast to the decrease in wetland area, CH₄ emissions from wetlands were similar at LGM and the preindustrial, with a value of ∼110 Tg yr⁻¹ (Figure 3). During the late Pleistocene, CH₄ emissions reached a maximum. At all time slices, however, global CH₄ emissions are within the range estimates for the present day [Kaplan, 2002], although seasonal and interannual variability may have been much more pronounced at times before the PIH (J. O. Kaplan, unpublished results, 2005). Global CH₄ emissions are dominated by tropical wetland sources in these simulations, and so the evolution and submergence of boreal wetlands during the course of the deglaciation only has a secondary effect on the global total.

[37] The late Pleistocene maximum in wetland CH₄ emissions occurs when a significant part of the tropical continental shelves were still exposed, while increasing precipitation due to warmer sea surface temperatures, and increased atmospheric CO₂ concentrations had a strong effect on vegetation productivity and carbon turnover. As sea level stabilized near present levels, and shallow continental shelves were submerged, CH₄ emissions from wetlands gradually decreased to near-preindustrial levels.

[38] Uncertainty in estimates of global wetland CH₄ emissions for the present day is large. While the total budget of CH₄ is well constrained by atmospheric measurement, a twofold range in the magnitude of global wetland CH₄ sources has been proposed by different studies [see Kaplan, 2002], though most studies converge in a range of ±20% of the mean. In the past we may expect somewhat greater uncertainty. While the processes of CH₄ formation in wetlands are identical, long-term changes in soils, topography, and CO₂ concentration may all effect net CH₄ emissions. We estimate that a doubling or halving of the global CH₄ emissions in our paleoscenarios would be within the range of uncertainty in model parameterizations and assumptions, although we do not expect the direction of change relative to PIH to be different. Future research should attempt full sensitivity analysis to better quantify the impact of these uncertainties on the composition of paleoatmospheres.

3.3. BVOC Emissions

[39] Table 2 summarizes BVOC surface flux simulated by BIOME4-TG at the four key time periods considered in the chemistry-climate model runs, along with estimates of BVOC emissions from other recent studies. In our simulations, total global BVOC emissions increased steadily from 21 ka to the preindustrial, with an LGM minimum of 625 Tg yr⁻¹ and a preindustrial value of 1050 Tg yr⁻¹ (Figure 3). BVOC emissions do show an early Holocene maximum with values at 10 ka near those for the preindustrial, a slight decrease to 6 ka, followed by recovery to near present-day levels.

[40] Our estimates of total global BVOC flux are consistent with Adams et al. [2001] and Valdes et al. [2005], though they differ by up to 30% for individual species. All estimates virtually agree on preindustrial emissions of isoprene and terpenes, but diverge more pronouncedly at the earlier periods, with the calculations of Adams et al. [2001] 15% to 30% higher, and Valdes et al. [2005] lower by similar magnitude. Larger differences arise for ORVOC emissions, which in our case are somewhat lower than the other estimates (between 25% at PI and 75% at 10 ka).

[41] Because of the short lifetime of BVOCs, the geographic distribution of emissions is as important as the total flux and the species composition. Over the last 21 ka the global pattern of BVOC emissions changes (Figure 5) for two major reasons: (1) the establishment of temperate and boreal vegetation after deglaciation and (2) the in situ enhancement of tropical sources stimulated by increased humidity and CO₂ fertilization.

[42] Ice sheet retreat and the subsequent development of cold needleleaf forests (Figure 2) resulted in the appearance of BVOC emissions from areas that had little or no emissions at the LGM. Temperate forest development also resulted in increased BVOC fluxes over eastern North America, Europe, and northeast Asia. Tropical BVOC emissions, which are significantly larger in overall magnitude than boreal or temperate sources, largely intensified in situ in response to the increased humidity and CO₂ concentrations that accompanied deglaciation. As noted in previous sections, the increase in some areas of tropical forests (e.g., in South America) was offset by the inundation of other regions (e.g., on the Sunda Shelf), which led to a meridional redistribution of tropical BVOC sources.

[43] The geographic changes in BVOC emissions are also apparent in the temporal pattern of the species distribution of total BVOC (Figure 6). Isoprene, which is largely emitted by broadleaf vegetation in tropical and temperate forests, increased most strongly between 21 ka and 10 ka, when the effect of CO₂ fertilization and increasing humidity would have been strongest. In contrast, terpene emissions, which are greatest in needleleaf forests such as those in cold biomes, increased steadily to the present because of the slow development of boreal forests in deglaciated regions, and somewhat less so by the nonuniform effect of warmer temperatures in the high latitudes. Methanol fluxes (included in other BVOC) originate primarily from grasslands [Galbally and Kirstine, 2002; Kirstine et al., 1998; Warneke et al., 1999], and also increased over the entire period of deglaciation and throughout the Holocene. Increase in these fluxes is related to the increase in the productivity of global grasslands as a result of climate warming and increased atmospheric CO₂ concentrations.

[44] In all of our calculations of BVOC emissions substantial uncertainties exist, not in the least because of the uncertainty in present-day estimates of VOC emissions rates from global terrestrial vegetation. Simulations of the global BVOC surface flux by vegetation models is affected by uncertainty in our simulation of vegetation distribution and climatic parameters controlling BVOC emissions (e.g., surface solar radiation and canopy deposition), as well as biological factors not considered in the vegetation model (e.g., herbivory and disease) [Guenther et al., 1995]. Because of the large measured plant-interspecific variation in BVOC emissions, the use of plant functional types in our model represents an additional simplification and adds to the uncertainty in our calculations. The assumption that VOC emission rates at conditions prevailing at the LGM and the interglacial periods are comparable to emission rates derived for present-day conditions remains unverifiable.

[45] For the present day, landscape-averaged emission factors for isoprene and monoterpene in similar ecosystems vary ∼fivefold [Guenther et al., 1995]. OVOC emission factors are much more uncertain, and could vary up to 100 times among species [Guenther et al., 1995; MacDonald and Fall, 1993]. This wide variability in measured emission factors could lead to uncertainty in the total global BVOC flux on the order of ±100%. On the other hand, plausible upper and lower limits for global BVOC emissions may be estimated from measurements of key atmospheric species, for example, O₃ [Folberth et al., 2005a]. Given relatively well-constrained budgets of well-mixed reactive trace gases for the present day [see, e.g., Ramaswamy et al., 2001], we believe that our calculations are a valid approximation of the relative changes in total annual BVOC emissions from global terrestrial vegetation from the LGM to preindustrial times.

3.4. OH Concentrations and CH₄ Lifetimes

[46] Table 3 summarizes the global annual mean burdens of CO, tropospheric O₃ and OH, and the CH₄ lifetime simulated by LMDz-INCA at each of the four model time slices. The decrease in OH concentration from the LGM to PIH of 2.8 × 10⁵ molecules/cm³ (22%) is concomitant with increases of similar magnitude in CO and O₃ burden. Because surface emissions of CO were kept constant in the three paleo time slices in these experiments, apparent changes in CO burden are related to the changes in OH. These results, along with the decrease in O₃ over the same time period, illustrate the effect of changing BVOC emissions on the oxidizing capacity of the atmosphere from the LGM through the Holocene.

[47] The vast majority of CH₄ and BVOCs emitted by the terrestrial biosphere are oxidized in the gas phase through the reaction with OH. An increasing reactive carbon burden in the atmosphere will reduce OH abundance, especially in the highly NO_x limited atmosphere of the preindustrial eras. Given a relatively constant source of CH₄ as simulated by BIOME4-TG, the change in the ratio of CH₄ to other RCCs drove the increase in CH₄ lifetime from LGM to preindustrial. This CH₄-to-other-RCC ratio decreased from 0.25 at LGM to about 0.19 at preindustrial. Correspondingly, the global annual mean photochemical lifetime of CH₄ increased by 2.1 years or about 29% from LGM (7.3 years) to preindustrial (9.4 years). The change in CH₄ lifetime is directly related to changes in its atmospheric concentration, as described in the following section.

[48] As noted above, our calculations of the changes in OH burden are based only on the change in surface emissions of reactive trace gases and do not include the effect of changes in climate, surface albedo, the stratospheric ozone column, or emissions of other reactive trace gases. These limitations are discussed below.

[49] Changes in the Earth's albedo from, for example, the presence or absence of sea ice and ice sheets, will have an effect on the global radiation budget. Ultraviolet (UV) radiation flux in the 290 to 310 nm range controls tropospheric ozone photolysis and therefore has a strong influence on OH burdens. However, studies by Pinto and Khalil [1991] and Martinerie et al. [1995, Figure 8] suggest that variations in the total ozone column between preindustrial and LGM were less than 2%. Even assuming that variations in stratospheric ozone were negligible [Pinto and Khalil, 1991], the impact of changing albedo on tropospheric chemistry is expected to be much smaller than the effect of changes in water vapor.

[50] The change in the water vapor content of the troposphere is a major contributor to the uncertainty of our model calculations. By repeating LGM model runs with preindustrial climate, Valdes et al. [2005] report a climate sensitivity of the global mean CH₄ concentration and OH burden on the order of 10% and 7%, respectively. Similarly, 2D chemistry-climate model experiments by Martinerie et al. [1995] suggested that lower water vapor during the LGM would decrease global OH concentration by 7%, resulting in an increase of CH₄ lifetime by 6%. On the other hand, uncertainty in regional climate change itself between LGM and PIH remains large, particularly for important parameters that influence tropospheric water vapor content, for example, tropical temperatures [Farrera et al., 1999; Pinot et al., 1999].

[51] Another source of uncertainty in the chemistry-climate model experiments is the potential for long-term changes in the production of NO_x by lightning. To investigate this possibility, we analyzed global fields of convective precipitation rate (a diagnostic for lightning abundance) in the output of the HC-UM paleoclimate simulations. Over the full 21 ka, no significant trend in convective precipitation rate could be identified; we therefore expect any changes in lightning-produced NO_x to have been small. Furthermore, the sensitivity analyses of Valdes et al. [2005] and those in previous studies [Sinha and Toumi, 1997; Thompson et al., 1993] indicate that variation in NO_x appears to be a minor component in influencing the OH and CH₄ budget.

[52] Finally, the direct effect of tropospheric O₃ concentrations on plant productivity, while well known [Gregg et al., 2003], was not included in this study. It is conceivable that reduced O₃, simulated here at all times before the PIH, could have led to an enhancement of overall vegetation productivity and hence BVOC emissions. While the effect is probably small, investigation of such feedback would be an important question for future research.

[53] Changes in climate and land surface conditions could counteract up to one third of the 22% change in OH burden between LGM and PIH we simulate here. Nevertheless, it is clear that changing BVOC emissions have a significant impact on tropospheric chemistry and CH₄ lifetimes. Further research and sensitivity analysis with 3D chemistry-climate models will be required to better quantify the importance of individual uncertainties, all of which contribute to the mismatch between modeled and measured CH₄ concentrations.

3.5. Modeled and Ice Core Measured Atmospheric CH₄ Concentrations

[54] In Figure 7 the atmospheric CH₄ concentrations simulated by LMDz-INCA are plotted together with measured CH₄ concentrations sampled from several ice cores from Greenland and Antarctica. Measured atmospheric CH₄ concentrations increased by about 450 ppb from LGM to preindustrial, with a rapid increase during the deglaciation, and a relative minimum in the mid-Holocene (∼6 ka). Modeled atmospheric CH₄ concentrations increase ∼385 ppb between LGM and PIH and are within 15% of the measurements at all time slices.

[55] We expect discrepancies between measured and modeled CH₄ concentrations in these experiments. As described in previous sections, each component of the linked model simulations that ultimately govern CH₄ concentrations has errors related to uncertainties and assumptions in the model setup. Nevertheless, the close agreement between model and measurements in these experiments demonstrates that changing BVOC emissions are a decisive factor in the CH₄ cycle that is, for example, more important than the direct effect of climate change on photochemistry and OH burden. Our results support recent hypotheses that Glacial-Interglacial changes in atmospheric CH₄ concentrations could have been driven primarily by changes in the photochemical sink for CH₄ rather than changes in the source [Adams et al., 2001; Valdes et al., 2005].

[56] The Holocene record of atmospheric CH₄ concentrations bears further examination. Despite the relatively stable climate over this time period, CH₄ concentrations show a marked long-term trend, with an early Holocene maximum followed by a minimum around 6 ka and a steady increase to PIH. Modeled CH₄ concentrations are higher at both 10 ka and PIH than 6 ka. The early Holocene maximum in CH₄ concentrations is a reflection of the maximum in wetland CH₄ emissions, which decreased ∼20% between 10 ka and 6 ka. At 6 ka, lower wetland CH₄ emissions effected lower atmospheric CH₄ concentrations despite the steady increase in BVOC emissions, as evidenced by the relatively stable OH burden over this time period (Table 3).

[57] Between 6 ka and PIH, the 15% and 30% respective anthropogenic increases in the CH₄ and CO inventories used in our experiments apparently have an effect on atmospheric CH₄. These results may support the recent hypothesis that early anthropogenic activity was a driver for the observed increase in atmospheric CH₄ concentrations over the past 6 ka [Ruddiman and Thomson, 2001]. It is likely, however, that areas of prehistoric rice cultivation were on floodplain wetlands, and that areas of anthropogenic CH₄ emissions supplanted areas of natural wetland emissions. This effect was not considered in this study and further research is needed to quantify the potential effect.

[58] Unfortunately, in these simulations we were not able to test mechanisms to explain rapid changes in atmospheric CH₄ observed in ice cores. We speculate that time-discrete mechanisms, such as the rapid inundation of tropical continental shelves and in situ decay of organic matter, concomitant with the long-term changes in the CH₄ sink, could explain the observation of rapid but permanent changes in atmospheric CH₄ concentrations. Further study with fully coupled land surface atmospheric chemistry models will be necessary in order to resolve remaining questions regarding the mismatch between model and observations.

3.6. CH₄ Interhemispheric Gradients

[59] Evaluation of the interhemispheric gradient in atmospheric CH₄ at times in the past is possible because of the existence of ice cores in the high latitudes of both hemispheres. These interhemispheric gradients may be used to imply temporal changes in the geographic distribution of CH₄ sources and sinks. In contrast to CH₄, which is relatively well mixed in the atmosphere, the short lifetimes of BVOCs mean that their effect on OH is mostly local. Therefore it is instructive to compare model simulations of the interhemispheric gradient in CH₄ to ice core records as a way of deducing the temporal variability in the interplay of CH₄ sources and sinks.

[60] Figure 8 shows normalized interhemispheric gradients of zonally averaged atmospheric CH₄ concentrations at each time slice of our model simulations and for the present-day simulation of LMDz-INCA (note the different scale in Figure 8 for the CH₄ gradient at present-day conditions because of the predominance of anthropogenic sources). The gradients were derived by calculating zonal means of the annually averaged CH₄ concentration at the lowest model level at each time slice and normalizing these values by subtracting the minimum CH₄ concentration for each time slice individually. In all cases, the minimum in atmospheric CH₄ occurs at approximately the same latitude (∼70°S). The appearance of the minima at this latitude is due to the absence of CH₄ sources in the high latitudes of the Southern Hemisphere. Poleward of the minimum, very low OH concentrations yield somewhat longer CH₄ lifetimes and hence marginally higher concentrations. Northward of the minimum, both surface fluxes and OH burdens increase, though with substantially different patterns at the different time slices.

[61] At LGM, the CH₄ gradient is nearly symmetric around the equator with a strong peak in the tropics due to the predominance of CH₄ sources in tropical latitudes (Figure 4). Low temperatures and widespread areas of ice, tundra, and cold forests (Figure 3) led to relatively low extratropical BVOC production, resulting in short CH₄ lifetimes and low concentrations at northern and southern midlatitudes to high latitudes. The interpolar difference (IPD) [Dällenbach et al., 2000] in CH₄ concentrations simulated by the model is 4 ± 4 ppb, which is comparable to the ice core LGM mean of −3 ± 4 measured by Dällenbach et al. [2000] and somewhat smaller than the 18 ± 7 reported by Brook et al. [2000], though it must be noted that both of these measurements represent integrations of some thousand years around 21 ka.

[62] During the late Pleistocene and early Holocene the CH₄ gradients show an increasingly asymmetric shape with an evolving second peak at northern midlatitudes to high latitudes (∼60°N). We attribute the evolution of the asymmetry in the CH₄ gradient over time to both the development of boreal wetlands after deglaciation, and the increase in BVOC emissions in the Northern Hemisphere and the consequent reduction in the strength of the photochemical CH₄ sink. Measurements of the IPD over this time period were approximately 44 ± 4 ppb [Chappellaz et al., 1997; Dällenbach et al., 2000], while the model simulates a substantially smaller IPD of only 14 ± 6 ppb. The origin of this model-data mismatch may be due to the omission of circulation-related climate changes in the LMDz-INCA model simulations and highlights a limitation of our modeling approach. The transport simulated by LMDz-INCA has been extensively validated for present-day conditions, both with passive tracers and for CH₄ itself [Folberth et al., 2005a; Hauglustaine et al., 2004]. However, different boundary conditions between PIH and 10 ka and 6 ka, particularly regarding solar forcing and the intrahemispheric temperature gradient, may have resulted in significantly increased latitudinal atmospheric transport relative to the PIH. While it is also possible that we underestimate the magnitude and development of boreal CH₄ sources in BIOME4-TG, recent syntheses of global methane sources indicate the general predominance of wetland CH₄ emissions from the tropics [Frankenberg et al., 2005].

[63] At the PIH, northern midlatitudinal and high-latitudinal CH₄ concentrations surpass the extant tropical peak in magnitude. Additional anthropogenic CH₄ sources from landfills and rice cultivation, along with a decreasing sink for CH₄ effected by increasing BVOC emissions, added to the buildup of the Northern Hemisphere midlatitudinal maximum.

[64] Under present-day conditions, the tropical maximum appears to disappear entirely because of the abundance of anthropogenic sources of CH₄ and photochemical CH₄ precursors, which are strongly biased to the Northern Hemisphere. The shape of the interhemispheric gradients in the Southern Hemisphere suggests that southern hemispheric CH₄ concentrations are dominated by atmospheric transport rather than direct sources.