Modelling holocene carbon accumulation and methane emissions of boreal wetlands – an earth system model approach

Introduction Conclusions References


Introduction
Wetlands, and particularly peatlands of the boreal latitudes, store considerable amounts of carbon (C) in the form of peat and constitute a significant natural source of methane (CH 4 ) even though they cover only 3 % of the global land surface. Previous studies suggest that the size of the boreal peat carbon stock is as large as 180 to 621 Petagramm C (Pg C) (Gorham, 1991;Smith et al., 2004;Turunen et al., 2002;Yu et al., 2010) and CH 4 emissions are in the range of 32 to 112 Teragramm per year (Tg CH 4 yr −1 ) (Bergamaschi et al., 2007;McGuire et al., 2010;Zhuang et al., 2004). By storing carbon and taking up carbon dioxide (CO 2 ) from the atmosphere, boreal peatlands have had a cooling effect on climate during the last millennium (Frolking and Roulet, 2007). This positive carbon balance also applies to the greenhouse gas (GHG) emissions, in particular CH 4 emissions of boreal wetlands. Undisturbed boreal peatlands are likely to continue functioning as a net carbon sink Tolonen and Turunen, 1996). On the other hand these carbon pools might be destabilized in the future since they are sensitive to climate change (Christensen, 1995;Dise, 2009;Kayranli et al., 2010). Given that the processes of peat accumulation and decay are strongly dependent on hydrology and temperature, this balance may change significantly in the future. Considering the projected future warming, boreal peatlands could potentially have a large impact on carbon cycle-climate feedback mechanisms and therefore play an important role in global carbon cycle dynamics (McGuire et al., 2009). However, global biogeochemistry models used for simulations of carbon cycle dynamics in past and future climates usually neglect peatland processes (Frolking et al., 2009).

Modelling carbon cycling in boreal wetlands
The boreal wetlands existing today were established after the Last Glacial Maximum and have continued to grow during the Holocene (Jones and Yu, 2010;Yu et al., 2010). They are linked to the terrestrial carbon cycle in many respects. On the one hand there is the large carbon stock in the soil which was built up despite the comparatively low net primary production (NPP) with an average of 100 to 400 g C m −2 yr −1 (Blodau, 2002). On the other hand these boreal wetlands produce carbon emissions, which, as a result of several competing processes, may exceed the carbon uptake and turn the peatland from a carbon sink into a carbon source. These emissions occur either in the form of CO 2 due to oxic decomposition or in the form of CH 4 due to anoxic decomposition of organic material. As CH 4 has the second-largest radiative forcing of the long living GHGs after CO 2 , it is of particular importance to identify the composition of carbon emissions from boreal wetlands.
As the key features of these northern ecosystems -carbon accumulation and GHG emission -will be investigated in this study, these ecosystems will be referred to hereinafter as boreal wetlands as the generic term, or peatlands, indicating the accumulation of carbon in the soil. Our term "boreal wetlands" only refers to permanent wetlands of the northern latitudes, but not to seasonally inundated areas. In the current literature, however, three main terms are used to encompass the subject: wetlands, peatlands and mires (Rydin et al., 2006).
In the soil column of peatlands there is a sharp transition between the upper and the lower soil layers (Charman, 2002), with different biogeochemical properties. The reason for this is the high water table (wt), which divides the soil column into an oxic part above and an anoxic part below the water table. Ingram (1977) named the upper layer acrotelm, and the lower one catotelm. Usually, the boundary between the acrotelm and the catotelm is defined as the deepest point to which the water table descends in an annual cycle (Charman, 2002).
Wetlands form peat when conditions are suitable for the growth of plants and hydrology inhibits the aerobic decomposition of dead organic matter. Therefore, organic matter accumulation in wetlands is a function of the balance between net primary productivity and abiotic and biotic decomposition processes (Reddy and DeLaune, 2008). The biogeochemistry of wetland soils is primarily controlled by hydrology, which in turn is first and foremost represented by the water table. Carbon that is fixed in the vegetation becomes litter when plants die, either from disturbances or at the end of the vegetation period. Later, in the acrotelm, the main plant structure collapses and becomes part of the catotelm, where it finally accumulates. These structural layers of lit-ter, acrotelm, and catotelm, in principle, are fixed and continuously move upward with the growth of the peat surface (Charman, 2002). The boundary of the functional zones between acrotelm and catotelm, however, does move, particularly when the position of the water table changes through the seasons. A key factor in determining chemical transformations in peatlands is the degree of aeration (Charman, 2002). The decomposition of the acrotelm happens relatively fast. The catotelm is anaerobic during the entire year; therefore, only anaerobic microbes decompose the organic matter very slowly.
Carbon emissions from boreal wetlands consist either of CO 2 or of CH 4 . Two main mechanisms govern the amount and ratio of CO 2 to CH 4 that is emitted to the atmosphere: again the degree of aeration, and the microbial CH 4 production and its oxidation (Kamal and Varma, 2008;Sundh et al., 1994). The production depends on the composition of the microbial community and several abiotic factors such as the availability of suitable organic material, soil temperature, and soil moisture. The quantity of CH 4 being oxidized in the soil depends strongly on the pathway of the freshly produced CH 4 to the surface (Dinsmore et al., 2009;Wille et al., 2008). Three pathways are common in boreal wetlands: diffusion, plant mediated transport, and ebullition. Diffusion through the soil column could lead to a strong oxidation of CH 4 whereas plant mediated transport through the stems of the plants reduces the probability of the oxidation of CH 4 . Ebullition, or bubble formation, leads to little oxidation since the bubbles rise quickly to the surface.
Apart from the functional segmentation in the vertical plane, peatlands are characterized by environmental gradients in the horizontal plane: patterns of vegetation variation can be found everywhere in peatlands. The change from floating mats of sphagnum or sedges at the edge of a water body, through taller vascular plants away from the water, and perhaps to a shrub community and a forested margin at the edge of a peatland, is obvious. Many of these patterns are brought about by gradual spatial variations in environmental conditions (Charman, 2002). The largest boreal wetland areas are located between 50 • N and 70 • N and are classified as bogs and fens (Reddy and DeLaune, 2008).

Local peat models
The first site-specific models of peat bog growth were developed in the late 70s and early 80s by Wildi (1978) and Clymo (1984). These site-specific models describe the growth of peat as a dynamic imbalance between input of plant material and its decay in the soil column. Clymo (1984) used two distinct carbon pools to represent the peat dynamics which constitute the functional layers of the acrotelm as the "active" zone and the catotelm for storing the carbon as defined by Ingram (1977). During the following years a number of peat models were developed with improved representation of peat accumulation processes. The geophysical model by Ingram (1982) takes into account the interaction of peat with the surrounding waters and illustrates how peat deposits are formed. Hilbert et al. (2000) analyzed the positive feedback between water table and peat accumulation in bogs in depth. They found two possible equilibrium states: one with deep water tables in drier sites, where peat depth increases with increasing water input, and another one for wet conditions, where the water table is near the surface. Annual litter cohorts were used as an input in the model by Bauer (2004) to investigate the effect of different litter quality on the peat accumulation. She found that a different vegetation community above a peatland could alter the peatland response to climate change. Yurova et al. (2007) used the general ecosystem model GUESS-ROMUL to investigate the effect of the water table depth on net CO 2 uptake rates of a sphagnumdominated mire. Their model simulations indicate that small changes in the water balance only have minor effects on the annual rates of land-atmosphere CO 2 exchange, whereas a water table higher than −10 cm or lower than −20 cm lead to a decreased net CO 2 uptake of the peatland. Based on an approach of modelling annual peat cohorts, the Millennia model by Heinemeyer et al. (2010) comprises climate-driven water table dynamics with a parameterization constrained by pollen-based vegetation reconstructions. The latest advance of local peat models is represented by the approach of Frolking et al. (2010), who described the accumulation and decay of peat with a definition of annual peat layers.

Global peat models
Recently, Kleinen et al. (2012) presented a study where both boreal peat growth and wetland extent are modelled in combination within the Dynamic Global Vegetation Model (DGVM) LPJ. They showed an accumulation of 330 Pg C during the last 8000 yr for areas above 40 • N and a modelled wetland area that scores well in comparison with data. This approach allows for modelling peatlands in different climate states, such as previous interglacials, since all required inputs are determined internally. Wetland fraction and water table from this approach were used in our model setup.

Methane emission models
The importance of peatlands with respect to GHG emissions was taken into account by the development of a processbased model by Cao et al. (1996), who drove their wetland methane emission model with data sets for climate, vegetation, soil, and wetland distribution. It could be shown that the amplitude of CH 4 emissions depends on a combination of available soil carbon, its decomposition rate, soil moisture and bacterial activity.
The first process-based CH 4 emission model that was calibrated for specific sites was developed by Walter et al. (1996). It was first published as a one-dimensional model, but has been developed further for global applications (Walter and Heimann, 2000;Walter et al., 2001a, b). It describes the production and oxidation of CH 4 in the soil column and accounts for three pathways of CH 4 to the surface: diffusion, plant-mediated transport, and ebullition. This mechanistic approach has been a starting point for other CH 4 models. Van Huissteden et al. (2006) implemented this approach into a site-specific model to assess the effect of water table management on CO 2 and CH 4 fluxes from peat soils. Wania (2007Wania ( , 2009aWania ( , b, 2010 implemented Walter's model in the DGVM LPJ. Furthermore, they included processes such as permafrost and wetland specific plant functional types to investigate the CH 4 cycle in boreal regions. The key benefit of this approach is that -except for climate data, a soil map, and prescribed wetland extent -all input data are determined by the DGVM itself. Gedney et al. (2004) developed a simple CH 4 emission scheme that runs within the UK Met Office climate model HadCM3. They parameterized the flux from wetlands by including the basic controls of temperature, water table position and soil carbon. Eliseev et al. (2008) implemented a module of CH 4 emissions from wetland ecosystems and a module for soil thermal physics into their climate model of intermediate complexity. The extent of wetlands was prescribed. They showed that CH 4 emissions increased from 130-140 Tg CH 4 yr −1 (pre-industrial) to 170-200 Tg CH 4 yr −1 at the end of the 21st century. Meng et al. (2012) and Riley et al. (2011) developed a biogeochemical model and integrated it into the Community Land Model (CLM4Me) with the purpose of understanding the uncertainty and its sources that emerged during the development and application of the models described above. Because of their advanced model of CH 4 production and emission, they found large sensitivities in CH 4 fluxes to changes in model parameters and express low confidence in the predictions of future terrestrial CH 4 feedback strength. For an overview of peat and CH 4 emission models, see Table 1.

Outline
This study aims at the evaluation of the effect of boreal wetlands on climate through fluxes of CO 2 and CH 4 on millennial timescales, based on plausible peat accumulation patterns (Sect. 1). We developed a generic model of peatland carbon dynamics embedded in an Earth system model. Given that our intention is to run this model globally over long timescales such as the Holocene, we pursue the strategy of using a simple model that captures the main processes, i.e. plants produce litter which successively enters the acrotelm and finally, if not respired, becomes a part of the catotelm. This is important for the land carbon balance throughout  Wildi (1978) Interaction between plants, Functions derived from .
x x x x water and peat field meas. Clymo (1984) Peat growths model with Productivity and decay x x x x proportional decay function rates Hilbert et al. (2000) Peatland dyn. Interactions Precip., wet and dry x x x x between water table and peat growths conditions Yu et al. (2001) Generic model framework Water budget data x x x x Frolking et al. (2001) Long-term peat accumulation Vegetation, NPP, decomp.
x x x x model rates Bauer (2004) Litter-quality dependent peat NPP, water table position x x x x dynamics over millennia Ise et al. (2008) Physical We compare the results of this model with observations and results from an inversion model for the present-day situation (Sect. 3). This allows us to address the question of changes in CH 4 emissions from boreal regions and their contribution to the atmospheric concentration reconstructed from the ice cores at millennial timescales (Sect. 4). We end with a conclusion in Sect. 5.

Model description
JSBACH is the modular land surface scheme of the MPI-ESM (Raddatz et al., 2007;Roeckner et al., 2003). Within JSBACH the CBALANCE model describes the changes in carbon storage from the growth and death of plants and the remineralization of carbon in soils. Since we focus on the biogeochemical cycles of boreal wetlands, we modified the CBALANCE submodel by integrating wetland specific parameterizations. The modified CBALANCE submodel is henceforth called peatBALANCE.
Carbon accumulation due to the slow decomposition of plant material under anoxic conditions is not represented in the CBALANCE model. To investigate the carbon cycle dynamics of boreal wetlands, we extended the CBALANCE model to the typical wetland processes which have an effect on the carbon cycle: (i) the division of the soil column into the two functional layers of acrotelm (oxic and anoxic conditions during the course of the year) and catotelm (permanently anoxic conditions), (ii) a moving water table which defines these functional layers, and (iii) as an outcome of the interaction of the components listed above, the accumulation of carbon in the catotelm. Carbon emissions in the form of CH 4 are considered in this model setup by integrating the Walter model (Fig. 1).

The standard CBALANCE submodel in JSBACH
The CBALANCE model uses three different pools (green, litter, soil) to describe the storage of organic carbon in living, dead, and decaying plants, which are the state variables of the model. Generally, the land biosphere grows by filling the vegetation carbon pool with the carbon gained from photosynthesis (NPP G , Eq. 1). This "green pool" (C G , Eq. 1) contains carbon of the green or living parts of plants (leaves, fine roots, sapwood). Through seasonal leaf shedding, carbon is transferred from the C G to the "litter pool" (C L , Eq. 2), which is described by the flux F litter (Eq. 5). Here the litter flux is determined by assuming that the "green pool" is loosely coupled to the leaf area index (LAI) (Eq. 4) and cannot exceed a maximum value of the green pool (C G ) that is proportional to the LAI, where sla is the specific leaf area. The C L pool loses carbon due to heterotrophic respiration  (R L , Eq. 6), calculated through a Q 10 model, where T is the temperature of the uppermost soil layer in • C (Eq. 6). This fraction of slowly decomposing parts of the plant (β L R L , Eq. 3) is transferred to the "soil pool" (C S , Eq. 3). Finally, depending on the turnover time of the "soil pool", part of the carbon is respired to the atmosphere (R S , Eq. 7), using the Q 10 model as well.

The peatBALANCE submodel in JSBACH
The peatBALANCE model uses the same carbon pools for vegetation, plant litter, and the associated carbon fluxes as described above (Eqs. 1 to 2 and 4 to 6). Newly introduced are an "acrotelm carbon pool" (C A , Eq. 8), with temporarily oxic and anoxic conditions, that replaces the "soil pool" and a "catotelm carbon pool" (C C , Eq. 9), with permanently anoxic conditions.
The key factor controlling the decomposition rates of the total peat column is the position of the water table. Decomposition under anoxic conditions is slower by more than an order of magnitude. The fraction of the "acrotelm pool" decaying under oxic or anoxic conditions (f ae , Eq. 11) is determined by the position of the dynamic water table in relation to the height of the acrotelm (h A ), applying acrotelm density ρ and fraction of carbon C f in the acrotelm (Eq. 10). The acrotelm with oxic conditions decays with a faster turnover time (τ A ae = 14.92 yr −1 ) than the anoxic part of the acrotelm (τ A an = 42.64 yr −1 ); the values used were derived from (Malmer and Wallén, 1993). An overview of all parameters used is listed in Table 2. The acrotelm loses carbon due to heterotrophic respiration (k A , Eq. 12); the remaining slowly decomposing carbon passes into the catotelm (β A C A , Eq. 9).
The "catotelm carbon pool" is the place where the carbon accumulates as soon as the peat accumulation flux β A C A is larger than the respiration flux (k A , Eq. 12). Similar to Clymo et al. (1998) a turnover time of τ C = 30 000 yr −1 was set for the catotelm. The anaerobic respiration k C is modelled through a Q 10 model (Eq. 13). Parameter Value Name/Units Description C X -mol (C) m −2 (canopy) Size of carbon pool X F X -mol (C) m −2 (canopy) s −1 Carbon flux to the carbon pool X k X -mol (C) m −2 (canopy) s −1 Decomposition flux from carbon pool X NPP X -mol (C) m −2 (canopy) s −1 Part of NPP allocated to carbon pool X R X -mol (C) m −2 (canopy) s −1 Respiration flux from carbon pool X β X -Beta Fraction of the respiration flux that goes into the carbon pool X γ X 1.7 -Factor relating leaf carbon to the carbon content of the whole pool X f ae -- All carbon respired by the oxic parts of the soil column enters the atmosphere as CO 2 . The carbon respired by the anoxic parts, however, is a potential CH 4 flux. Carbon respired under anoxic conditions can either exist in the form of CO 2 or CH 4 . Potential CH 4 emissions are scaled by the ratio of CO 2 to CH 4 emissions (f m = 0.4) observed empirically (Scanlon and Moore, 2000) and this potential flux (P CH 4 , Eq. 14) is transferred to the methane emission model. The peatBALANCE model has a daily time step.

Methane transport model
In order to simulate methane emissions we employ the widely used (e.g. Bohn et al., 2007;Ringeval et al., 2010;van Huissteden et al., 2006;Wania et al., 2009a, b) Walter model (Walter and Heimann, 2000;Walter et al., 1996Walter et al., , 2001b. It is a process-based model which explicitly simulates the three most dominant pathways of CH 4 to the surface: i. Diffusion. The molecular diffusion of CH 4 through the water-logged soil column is slow and occurs everywhere. The calculation of this diffusion flux is based on Fick's first law.
ii. Ebullition. If a certain concentration of dissolved CH 4 is exceeded, CH 4 bubbles form. These are trapped in the peat and when a certain pressure threshold is reached, ebullition occurs. This release of CH 4 to the atmosphere happens so rapidly that only small amounts of CH 4 are oxidized.
iii. Plant-mediated transport. The transport of CH 4 through aerenchymous plant tissue from the place of origin directly to the atmosphere is defined as plant-mediated transport. Bypassing the aerobic zone in the peat column is a very effective transport mechanism leading to little oxidation of CH 4 . The plant-mediated transport to the surface depends on the distribution of roots in the acrotelm and catotelm and on plant phenology.
The production of CH 4 is calculated within the peatBALANCE model and distributed to the root zone of the methane model, which has a vertical resolution of 1 cm. Oxidation of CH 4 is only possible in the aerobic part of the soil column above the water table.
We coupled this model asynchronously to the peatBALANCE model described above and performed time slice model experiments with both models combined. The CH 4 emission model is called once per model day.
Since the purpose of the CH 4 emission model is solely to estimate the relation between CO 2 and CH 4 emissions, there is no feedback to the peatBALANCE model.
The model output from the peatBALANCE model consists of carbon accumulation rates, peat height, fraction of  carbon above or below water table, respiration of the oxic and anoxic parts of the soil, and the amount of carbon stored in the acrotelm and the catotelm. The coupled methane emission model gives the ratio and amount of CH 4 to CO 2 emissions.

Experiments
The peatBALANCE model is driven by soil temperature of the uppermost soil layer (topSoilTemp) (Fig. 2, left), net primary production (NPP) (Fig. 2, middle), leaf area index (LAI) and soil moisture (α). These driving data were extracted from a transient simulation with the MPI-ESM covering the last 6000 yr with orbital forcing on a yearly basis following Bretagnon and Francou (1988) and greenhouse gas concentrations set to pre-industrial values (Fischer and Jungclaus, 2011). In a previous publication (Kleinen et al., 2012) we described a scheme to dynamically determine the peatland extent and water table, based on the TOPMODEL approach (Beven and Kirkby, 1979), implemented in the CLIMBER2-LPJ model. We determined the water table distribution within a grid cell from the grid cell mean water table and topographic information. The summer mean grid cell fraction with a water table at or above the surface is considered an area wet enough for peat to accumulate. For the present experiments, the grid cell peatland fractions (Fig. 3) as well as the position of the water table within the peatland fraction are prescribed (Fig. 4), as determined by the CLIMBER2-LPJ model in a pre-industrial control run. Kleinen et al. (2012) showed that changes in the peatland area are rather small over the time period 8000 years BP to pre-industrial, so that we neglect changes in peatland area in the present study.
Using these boundary conditions we conducted Holocene simulations (6000 yr PB to pre-industrial) of    In order to get reliable starting conditions of the soil carbon pools (excluding catotelm peat, which never reaches    We focused our analysis of peat accumulation and CH 4 emissions on the boreal zone above 40 • N, subdivided into: Europe (EU), Asia (AS), and North America (NA). For comparison with data we defined smaller wetland areas as defined in Winderlich et al. (2011) for the West Siberian Lowlands (WSL) and Pickett-Heaps et al. (2011) for the Hudson Bay Lowlands (HBL), see Table 3 for coordinates.

Results for peat accumulation
Results from the transient 6000 yr Holocene experiment show an accumulation of 240 Pg C in the permanently anoxic catotelm pool (Fig. 5, left). The investigated boreal wetlands cover a region almost circumpolar at latitudes from 40 • N up to 80 • N. The peatBALANCE model simulates the largest increase in peat carbon storage in the areas around the HBL, in Europe, and in the WSL (Fig. 5, right). These are the major peatland areas observed today.
We compare modelled catotelm peat accumulation rates against observations of catotelm peat accumulation rates compiled from Gorham et al. (2003), Kremenetzki et al. (2003), Beilman et al. (2009), andYu et al. (2010), as     . 6. Modelled peat accumulation rates for pre-industrial times and at different sites compared to a data set of site scale peat accumulation rate measurements compiled by Kleinen et al. (2012) with data from Gorham et al. (2003), Kremenetzki (2003, Yu et al. (2010). described in Kleinen et al. (2012) (Fig. 6). Measured carbon accumulation rates are slightly larger, but the general patterns of peat accumulation are represented well by our model. This discrepancy could possibly be explained by a measurement bias since often the deepest (oldest) part of the peatland is sampled (Korhola et al., 2010) and the model represents a substantially larger area than a local peat site.

Pre-industrial times
Several hotspots of CH 4 emissions can be detected at 6000 yr BP, 3000 yr BP and pre-industrial (Fig. 7): the Hudson Bay Lowlands and Newfoundland, Eastern Europe and Finland, the West Siberian Lowlands and the Far East (Manchuria). We use the last of our 30-yr time slice model run to calculate the global pre-industrial CH 4 emissions.
To evaluate model performance, we compare the CH 4 emissions against data sets from two specific regions: the West Siberian Lowlands and the Hudson Bay Lowlands. These two data sets integrate CH 4 emissions over a larger area, which allows direct comparison with our global model output. Site scale CH 4 flux measurements are likely biased since they are measured where emissions occur, whereas our model averages over a comparatively large area. To upscale methane emission rates from individual study sites to large study regions, extensive field knowledge of individual land cover classes in an investigation area is indispensable (Schneider et al., 2009), but not available for the whole boreal zone.
Our model results show CH 4 emissions for the West Siberian Lowlands of 6 Tg CH 4 yr −1 (Fig. 8, blue). Based on the atmospheric concentration data from four tower stations in this region (ZOTTO; Sasakawa et al., 2012;Winderlich et al., 2010), an atmospheric transport inversion results in emissions of 6.89 Tg CH 4 for the year 2009 (Winderlich, 2011). The seasonal cycle of our model compares well to the results of the data-driven inversion (Fig. 9). However, the modelled summer emissions are smaller (13.2 instead of 17.5 g CH 4 m −2 yr −1 in July) than the atmospheric data implies. The strongly simplified freeze-thaw processes in our model could explain the enhanced winter emissions. Another modelling study by Krupchatnikow and Krylova (2004) in-   (Fig. 8, red). Pickett-Heaps et al. (2011) published observational data from two stations, quantified by using the GEOS-Chem chemical transport model. Their best estimate for this area is 2.3 Tg CH 4 in the year 2008.

Holocene time slices
Apart from the pre-industrial simulation, we also analyzed the CH 4 emissions for another six time slices, going back in time until 6000 yr BP. Our model results indicate that the spatial patterns of CH 4 emissions (Fig. 2, right and Fig. 7) change in pre-industrial times compared to 6000 yr BP. We find a decrease in CH 4 emissions in Newfoundland and the high Eurasian Arctic, and an increase in the Hudson Bay Lowlands area.
Adding up the emissions for circumpolar boreal regions (above 40 • N) results in 49.3 (± 2.3) Tg CH 4 yr −1 at 6000 yr BP, which increases to 51.5 (± 2.75) Tg CH 4 yr −1 in pre-industrial times. In general, the boreal emissions show an increase over the last 6000 yr, with some variability (Fig. 10,  right).
We analyzed the development of CH 4 emissions in three sub-boreal regions to get a detailed picture of where the carbon emissions originate: (i) North America, (ii) Europe, and (iii) Asia, (Fig. 10, left; Table 3). (i) In the North American domain we find a small increase of the emissions around 3000 yr BP, which decline to the initial values in preindustrial times. (ii) In the Asian domain we find an increase of CH 4 emissions from 6000 yr BP (10 Tg CH 4 yr −1 ) to preindustrial times (11 Tg CH 4 yr −1 ). (iii) In Europe there is an increase in CH 4 emissions between 6000 and 4000 yr BP, followed by a stable phase until pre-industrial times.

Peat accumulation
Comparing the areas where data of soil organic carbon content suggest present-day peatlands (e.g. Tarnocai et al., 2009) to the wetland area determined by the model, we see that the modelled wetland distribution contains most areas where Tarnocai et al. (2009) show high carbon concentrations (Fig. 11). Since climate varies with latitude as a first approximation, we assume that the good fit to the latitudinal distribution also leads to a reasonable representation of peat accumulation, though the exact distribution and area of peatlands may not be represented by the model (Kleinen et al., 2012). The modelled peatland distribution, however, has the advantage of providing a consistent wetland distribution at the global scale, whereas global data sets based on measurements come along with different rates of accuracy and uncertainty. Many investigations address the size of the boreal wetland peat carbon stocks. One common approach is to determine the basal peat age and to measure the height of the peat to derive carbon accumulation rates Yu et al., 2010). Estimates for peat accumulation rates for the Holocene, which are derived from radiocarbon analysis of peat cores, have been published by Yu et al. (2009). The timeweighted average rate in their study is 18.6 g C m −2 yr −1 during the Holocene, based on 33 peat cores from northern peatlands.
Our model shows a slightly lower number for the average peat accumulation rate, which is 16.6 g C m −2 yr −1 . Overall, the changes of peat accumulation rates throughout the Holocene are rather small. This is in line with other Holocene peat accumulation studies such as Yu et al. (2009). Borren et al. (2004) found the peat accumulation rate to be as Fig. 11. Distribution of soil organic carbon contents (0-100 cm depth) based on the NCSCD published by Tarnocai et al. (2009), displayed in the model resolution.
large as 30 to 50 g C m −2 yr −1 during the last 6000 yr (even with a peak of up to a 100 g C m −2 yr −1 shortly before the year 4000 PB). Since their study concerns the southern taiga of western Siberia, higher values than those found in a global approach including the higher latitudes can be expected.
When comparing measured to modelled peat accumulation rates, it should be kept in mind that our model calculates the average for the whole wetland fraction of the grid cell, whereas the site studies mostly are located in the centre of the wetland, which could possibly lead to higher values.
The peatBALANCE model accumulates 240 Pg C within 6000 yr of simulation. Considering the initialisation and expansion of peatlands in the early Holocene between 11 000 and 8000 yr BP (Jones and Yu, 2010;Korhola et al., 2010;MacDonald et al., 2006;Yu et al., 2010) and an estimated size of the accumulated carbon from the last glacial until preindustrial times of 180 to 621 Pg C (Gorham, 1991;Turunen et al., 2002;Yu et al., 2010), the accumulation of 240 Pg C over 6000 yr in our simulation fits well into this pattern of peat accumulation. The estimated carbon stock, however, does not include permafrost carbon stocks in the northern cryosphere region, which store an amount of 1400 to 1850 Pg C (McGuire et al., 2010;Tarnocai et al., 2009).
Our results demonstrate that the temporal behaviour of carbon accumulation is fundamentally different for the acrotelm (the periodically aerobic upper part of the peat profile) and the permanently anoxic catotelm (below the acrotelm). The acrotelm carbon pool reaches an equilibrium state with a size of 16 Pg C after a comparatively short period of 150 to 250 yr. This roughly corresponds to a peat layer depth of 0.4 m, which is in line with estimates in the literature (e.g. Charman, 2002).
Under Holocene climatic conditions, the amount of the carbon in the catotelm, in contrast to the acrotelm, grows almost linearly over the 6000 yr of transient simulation without a significant indication of saturation. In this time frame, peatlands are therefore a non-equilibrium system with regard to the carbon cycle.
If productivity at the surface and the relative decomposition rates in the acrotelm and catotelm remain constant over a long period of time, the total amount of respired (mineralized) carbon increases with time. Ultimately, the total peat decay reaches a level equivalent to the production at the surface and there will be no further net accumulation of peat mass. This is the principle upon which Clymo's (Clymo, 1984) idea of a limit to peat growth is based (Charman, 2002). The net peat accumulation becomes more variable if changes in water table conditions are incorporated (Hilbert et al., 2000). In the Hilbert model, peatlands are capable of switching from net sinks to sources of carbon quite rapidly (Charman, 2002). In our model and for our time frame of interest, the Holocene, the terrestrial carbon cycle in peatlands is a non-equilibrium system. Presumably, it will take much longer than 60 000 yr for these peatlands to come to an equilibrium at which catotelm decay of old carbon balances additions of carbon from the intermittently aerobic acrotelm above, as in the Clymo model (Clymo, 1984;Gorham et al., 2003).
The evolution of peat carbon stocks under future warming scenarios is a potential area of application for this modelling approach. Generally, an increase in soil temperature will increase the rate of organic matter decomposition; conversely, when soil is saturated with water, anaerobic conditions slow down decomposition rates.

Methane emissions
The intensity of CH 4 emissions from natural wetlands are very uncertain because these emissions vary considerably in time and space (Frankenberg et al., 2005) Our results show emissions in the range of 49.3 to 51.5 Tg CH 4 yr −1 for boreal wetlands north of 40 • over the course of our Holocene experiment from 6000 yr BP to preindustrial. We identified the increasing carbon stock as the main driver of rising CH 4 emissions. Lower soil temperatures as existing in the forcing would lead to a natural decrease in emissions, which we do not see in our results.
Measurements along the Dome C Antarctic ice core and the GRIP ice core in central Greenland show a minimum in CH 4 concentration around the year 5000 BP followed by an increase until the present day (Blunier et al., 1995;Flückiger et al., 2002). Results from our model experiment indicate that CH 4 emissions from natural boreal wetlands contribute to this trend of rising atmospheric CH 4 concentrations as an increase of 1 ppb CH 4 yr −1 would correspond to a CH 4 increase in the global atmosphere of 2.83 Tg CH 4 yr −1 . However, the slowly rising natural emissions cannot explain the rapid CH 4 increase towards the end of the Holocene as indicated by ice core measurements.
Not all three examined subdomains show this increase. In the North American domain, the emissions decline towards pre-industrial and in Europe there is no clear trend visible in the mean values. Only the Asian domain shows a clear rising trend (from 10 to 11 Tg CH 4 yr −1 ) and the maximum values of the European domain rise from 17 to 19 Tg CH 4 yr −1 .
For two selected regions, the Hudson Bay Lowlands and the West Siberian Lowlands -two major areas of wetlands today -the modelled pre-industrial CH 4 emissions agree well with observations. Both model and data show the highest emissions in summer (JJA) and the lowest emissions in winter (DJF). The modelled CH 4 summer emissions reproduce well the data from Winderlich et al. (2011), but its seasonal shoulders are too broad. The modelled winter emissions are higher than observed because the freezing of the soil and the associated discontinuity of the methane transport are not represented by our model.
Modelling carbon cycling in the HBL area involves additional uncertainties. The HBL were basically shaped by the Laurentide Ice Sheet, which disintegrated about 7800-8000 yr ago, and have been rising isostatically since then (Glaser et al., 2004). These low-lying wetlands started to accumulate peat shortly after their emergence from the sea. Yu et al. (2010) indicate basal dates in the HBL area younger than 8000 yr BP. Our fixed wetland distribution does not account for this evolution of wetlands. Regarding our time frame of interest -starting 6000 yr BP -this evolution, however, is negligible.
Typically, the HBL are assumed to contribute 10 % to boreal wetland emissions (Pickett-Heaps et al., 2011). Since our model results show CH 4 emissions of ∼ 50 Tg CH 4 yr −1 for the boreal zone, we would expect emissions of about 5 Tg CH 4 yr −1 for the HBL accordingly. If we choose the HBL exactly as defined by Picket Heaps et al. (2011) we get emissions of 1.6 Tg CH 4 yr −1 . Since our coarse model does not resolve the HBL in detail, we would like to point out that our model produces methane south of the defined HBL area as well, which should be included in the HBL budget (Fig. 7).
Spaceborne CH 4 concentration measurements with the SCIAMACHY instrument on board the ENVISAT satellite show high CH 4 concentrations in the area of Manchuria, China, in 2003(Frankenberg et al., 2005, although not every year. Our model shows CH 4 emissions in Manchuria as well, but they were not investigated in more detail because no measurements exist.

Limitations of the model approach
There are some limitations of this modelling approach. First, the model has a coarse resolution. Second, the model is simplified. It is exactly these two restrictions, however, that allow modelling over long time frames and the integration in global models.
Some processes, which are observed in reality, are missing in the current version of our peatBALANCE model: sulphur deposition, permafrost, and sphagnum-associated methane oxidation. Boreal wetlands have always been strong CH 4 sources, accounting for 3 % to 9 % of the net land source of 552 Tg CH 4 yr −1 estimated by AR4 (Denman et al., 2007). With increases in CO 2 and temperature, and the associated increases in wetland productivity, CH 4 fluxes would be expected to increase. However, SO 2− 4 deposition (from industrial combustion) has the potential to divert substrate flow away from methanogens and thereby inhibit CH 4 flux to the atmosphere (Schimel, 2004). Permafrost processes in the northern boreal zone inhibit soil decomposition leading to a huge build-up of frozen organic carbon and hinders methane production. And last but not least, sphagnumassociated methane oxidation occurs ubiquitously across the globe (Parmentier et al., 2011;van Winden et al., 2010) and could reduce the amount of CH 4 significantly.
Modelling carbon cycling of boreal wetlands in the context of global Earth system models remains challenging. Wetlands are highly heterogeneous ecosystems both in spatial extent and at the process level. The main focus of future research could be a better representation of this diversity in the model, even if it is supposed to run at the global scale. The identification of the major processes and finding the simplicity in complexity are prerequisites for modelling the carbon cycle of boreal wetlands within a global climate model. For further progress in this direction, it is important to include a dynamic wetland model in JSBACH such as developed by Stacke and Hagemann (2012).

Conclusions
We developed the peatBALANCE model with peat accumulation and decay and implemented it in the land surface model JSBACH of the MPI-ESM. In a transient Holocene experiment, the model was driven by soil temperature, LAI, soil moisture, and NPP from MPI-ESM simulations by Fischer and Jungclaus (2011) and pre-industrial wetland fraction and water table levels from CLIMBER-LPJ simulations by Kleinen et al. (2012). The peatBALANCE simulation yields an accumulation of 240 Pg C in the boreal wetlands over the last 6000 yr. This carbon accumulation is an important long-term component of the carbon cycle and it is essential that it be accounted for in the simulations of the carbon balance on a millennial timescale.
We furthermore coupled the Walter methane emission model to our peatBALANCE model. Simulated preindustrial CH 4 emissions are 51.5 Tg CH 4 yr −1 for boreal wetlands above 40 • N, which is in the range of estimates of 15-112 Tg CH 4 yr −1 for boreal wetlands. The model run suggests that CH 4 emissions were slightly lower at 6000 yr BP than at pre-industrial. Emissions rise until preindustrial times with some variability. These dynamics in boreal CH 4 emissions resemble the minimum in the atmospheric CH 4 concentration around 5000-4000 yr BP with the following increase. The rising trend in CH 4 concentration over the last several thousand years, therefore, may be explained not only by anthropogenic factors (CH 4 emissions from landuse; Ruddiman, 2003), but also by natural processes.
Our approach does not only add processes in boreal wetlands that were important for carbon balance in the past, but also provides a framework for accounting for wetland response to future climate change. This is important, as on long timescales boreal wetlands are a significant component of the carbon and methane cycles that could either amplify or dampen human-induced global warming.