Soil temperature response to 21 st century global warming : the role of and some implications for peat carbon in thawing permafrost soils in North America

Introduction Conclusions References


Introduction
Northern peatlands cover ∼10 % of the land north of 45 • N and contain ∼500 Pg carbon (C), almost entirely as peat (Yu et al., 2010). This relatively large carbon pool within Correspondence to: D. Wisser (dominik.wisser@unh.edu @) one to a few meters of the land surface confers the important role of peatlands in the earth's carbon fluxes (Gorham, 1991). Changing environmental conditions could potentially impact the exchanges of carbon dioxide (CO 2 , Limpens et al., 2008), methane (CH 4 , Christensen et al., 2004), dissolved organic carbon (DOC, Frey and Smith, 2005) and nitrous oxide (N 2 O, Elberling et al., 2010) between peatlands and the atmosphere through a number of direct and indirect effects, and thereby change the role of peatlands in the global climate system. In particular, climate warming can directly affect C fluxes from peatlands. Ecosystem respiration has been found to increase in response to a 1 • C soil warming in a subarctic bog (Dorrepaal et al., 2009) and to vary depending on seasonal temperature patterns (Lund et al., 2010). Methane emissions are also strongly dependent on soil temperature (Treat et al., 2007;Kettridge and Baird, 2008). 122 D. Wisser et al.: Soil temperature response to 21st century global warming in the next decades (Stendel and Christensen, 2002; Arctic Climate Impact Assessment (ACIA), 2004; Marchenko et al., 2008) have raised questions on the potential release of carbon from peatlands as a result of global warming (Christensen et al., 2004;Turetsky et al., 2007). Anisimov (2007) predicted an overall increase in the seasonal thaw depth of 30 % by 2080, which resulted in a 20-30 % increase in CH 4 emissions from Russian peatland regions.
Soil temperatures are not only influenced by air temperatures but also by the presence or absence of an insulating layer of snow on the surface (Sturm et al., 1997). The response of soil temperatures to changes in air temperature strongly depends on the temporal distribution, depth, and density of the snowpack during the winter (Lawrence and Slater, 2010). The snowpack generally insulates the soil from the cold air and changes in the snowpack can therefore lead to more soil cooling (if the snowpack decreases) or warmer winter soil temperatures and permafrost degradation with a deeper snowpack (Christensen et al., 2004;Payette et al., 2004;Stieglitz et al., 2003).
Finally, it is the thermal and hydraulic properties of the peat itself that lead to different responses to air temperature in peat soil when compared to mineral soils that have prompted their explicit parameterization in land surface models (e.g. Beringer et al., 2001). Thermal properties of peat soil are strongly related to the peat soil moisture (Farouki, 1981;O'Donnell et al., 2009), which depends on the position of the water table, both of which are likely to change under future climate conditions. The insulating properties of peat are expected to preserve permafrost in peatland areas from severe degradation (Shur and Jorgenson, 2007;Yi et al., 2007). Currently, relict permafrost formed during the Little Ice Age is found in isolated patches within Canadian bogs despite warmer air temperatures (Halsey et al., 1995;Vitt et al., 1994). Additionally, peat depth is an important control on the active layer depth (Fukui et al., 2008;Yi et al., 2006). It is therefore imperative to consider the spatial distribution of those properties in attempts to model the future distribution of permafrost areas and the role of peatlands in the global climate system. A layer of organic material (<0.3 m) on top of mineral soil is a dominant feature of northern forest and tundra soils that are not classified as peatlands . This layer plays an important role in exchanges of water, energy and carbon and the representation of this thin organic layer in ecosystem models is therefore also important (Yi et al., 2009).
A number of attempts have been made to incorporate peatlands and organic soil layers in Earth system models (e.g. Rawlins et al., 2003;Wania et al., 2009) and these approaches have assumed a spatially uniform depth of the organic layer. Here, we present the first attempt to estimate soil temperatures in peatlands using a grid cell based framework with partitioning of the grid cell into a peat and non-peat fraction and explicitly taking into account the spatially varying depth of peat over a large spatial domain. This allows for an assessment of the impact of peat depth on soil temperatures and the relative differences in soil temperatures for peat soils and mineral soils and an assessment of the frozen and unfrozen peat volumes under changing climate. We use a thermodynamic model and soil moisture from a hydrological model to predict soil temperature in peat and mineral soil and estimate areas and thawed volumes of these soil fractions under future climate conditions. Assuming a carbon concentration and carbon content for peat and mineral soils, the total soil carbon in Northern North America that will change from being permanently frozen to seasonally thawed can be assessed.
We will first present an overview of the model and the data sets and then discuss the importance of peat depth for the prediction of peat soil temperatures by analyzing results from transient simulations of soil temperatures in peatlands in North America under predicted climate conditions for the period 2001-2100.

The GIPL 2.0 model
We used the modified Geophysical Institute Permafrost Lab model (GIPL 2.0;Marchenko et al., 2008) to simulate soil temperatures and freezing and thawing depths in peatland and non-peatland soils. Details of the numerical algorithm are provided elsewhere (Marchenko, 2001;Marchenko et al., 2008) and only a short overview is given here. GIPL calculates soil temperatures at different depths by applying a finite difference scheme to solve the 1-D heat transfer equation with phase change numerically and is therefore applicable to transient simulations. The version applied here computes soil temperatures at a 30 min × 30 min spatial resolution. GIPL takes into account spatially varying thermal properties of the soil layers and estimates the amount of unfrozen water based on those properties. The depth of the soil column and its importance for longer timescale models of soil temperature has been recognized (e.g. Alexeev et al., 2007;. The computational domain is therefore extended to a depth of 100 m at which a spatially uniform geothermal heat flux is set as a boundary condition. The upper boundary condition at the soil or peat surface is determined by the air temperature, or, if a snowpack is present, by a heat flux through the snowpack that is calculated taking into account the depth and density of the snow (Fig. 1). We applied the model with a time step of 24 h and computed soil temperatures at 114 computation nodes in the soil column, with space between nodes increasing with depth. The entire model domain (Fig. 2) covered a total area of 11.59 million km 2 in 8078 grid cells. In areas with peatlands, soil temperatures are calculated for peatlands and mineral soils separately. There was no lateral transport of heat between the peat and the non-peat fraction of a grid cell and no lateral transport of heat between neighboring grid cells. For both peat and non-peat calculations, we set the initial soil temperature, T ini [ • C], at depth z [m] as T ini (z) = T m G 0 z where T m is the mean annual air temperature [ • C] and G 0 [ • C m −1 ] is the geothermal gradient (set to 0.015). We computed time series of soil temperature for the period 2001-2100 after a 100-year spin-up period using contemporary (2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010) climatologies of temperature and precipitation.
The latent heat of the freeze-thaw phase change is implemented through the apparent volumetric heat capacity over a small temperature range near the freezing temperature of water (T f = 0 • C), C(z, T ) + δ(T )Q(z) at the freezing front T (z, τ ) = T f , where Q(z) is the latent heat of ice fusion or specific enthalpy of fusion (Q f = 334 [M J m −3 ]) and δ(T ) is the Dirac delta function having the following properties: where is the width of a half-interval ofC (volumetric heat capacity) andλ (thermal conductivity) smoothing range, set to 0.08 • C. The latent heat of ice formation-fusion is therefore represented through apparent volumetric heat capacity [δ(T )Q(z)] at the phase freezing front so that the heat transfer problem can be solved without explicit representation of freezing front coordinates. The smoothed volumetric heat capacityC [M J m −3 K −1 ] at temperature T [ • C] is approximated taking into account the volumetric heat capacity C vol of the solid materials (peat or mineral soil; see below), the soil water/ice content θ, the unfrozen water content θ u , and the latent heat of ice fusion Q(z): where C w and C ice are the heat capacity values for water and ice respectively (4.18 and 1.9 M J m −3 K −1 ). Unfrozen water θ u is the amount of water in the soil (θ) that remains unfrozen in the form of a thin layer around the soil particles. It can improve the transfer of heat in frozen soils (Farauki, 1981) and can retard the thermal response of the soil to air temperature changes (Romanovsky and Osterkamp, 2000). Taking into account frozen water in soil temperature models can therefore significantly improve the simulations. The amount of unfrozen water depends on the soil type and shows significant variations in space. We calculated unfrozen water for temperatures below a predefined threshold c (set to −0.01 • C), depending on the soil type as (Tice et al., 1976): where a and b are soil type specific shape parameters that were set to 0.1 and −0.7 for sand and to 0.5 and −0.7 for silt and clay soils. Unfrozen water for the mineral soil in a grid cell was computed as the weighted average of the unfrozen water for sand, silt, and clay. Unfrozen water content is typically very small in peatland soils (Romanovsky and Osterkamp, 2000;Kujala et al., 2008) and was neglected in our model.

Snow dynamics
Heat transport in snow is largely dependent on the microstructure of the snow (pore and grain distribution and size) which is difficult to quantify (Sturm et al., 1997), particularly for continental-scale simulations. The thermal conductivity of snow is therefore often calculated using empirical functions that are based on measurable properties of snow. We used the empirical equation of Proskuryakov (1999) (e.g. Ershov, 2004), which calculates the heat flux through the snow pack into the soil as a function of snow depth D s [m] and snow density ρ s [kg m −3 ]: where α 0 [set to 20.14 W m −2 K −1 ]is the averaged factor of convection heat exchange at the surface of snow and α is the snow thermal resistivity. The thermal conductivity increases with increasing snow density and is therefore higher later in the snow season.
To simulate spatially varying fields of snow density and snow depth, we implemented a simple snow accumulation/ snow melt module of the WBM plus hydrological model   (Wisser et al., 2010). In this model, it is assumed that precipitation falls as snow if the daily air temperature is below a temperature threshold [−1 • C] and accumulates a snow pack [mm snow water equivalent, SWE].
Melting from the snow pack occurs when daily air temperature is above 1 • C as a function of daily air temperature and rainfall (Willmott et al., 1985): where M s [mm] is the snowmelt, T m [ • C] is the mean daily air temperature, and P [mm] is the daily precipitation. The density of snow can change by an order of magnitude from wild, fresh snow to compacted, settled snow later in the season as a result of metamorphism. To calculate the snow depth during the entire snow season for each grid cell we assumed an initial snow density ρ s [kg m −3 ] of 150 at the onset of the winter (defined as the middle of the first month where the mean air temperature is below −1 • C) and a linear increase in snow density of 3 kg m −3 day −1 throughout the winter based on regional estimates of snow depth by Gray and Prowse (1993) so the depth of snow for each time step can be calculated as D s = SWE/ρ s .

Soil thermal properties
Soil apparent volumetric heat capacity C vol [J m −3 K −1 ] and soil thermal conductivity λ [W m −1 K −1 ]) are calculated as the volume weighted average values for the fractions of solid material (peat or mineral soil), and water. The depth of the organic layer on top of mineral soil depends on soil type, vegetation, and land cover. We used the 0.5 • version of the MODIS land cover map (Friedl et al., 2002) to assign organic layer depths to land cover classes ( Table 1). The thermal properties of the organic layer were assumed to be the same as the peat thermal properties. The total volume of this non-peat organic layer overlaying mineral soils in the model domain is 1173 km 3 , and its depth ranges from 0 to 0.25 and is 0.12 m on average.

Thermal properties of mineral soils
Heat capacity c for dry mineral soils was calculated using an empirical relationship between bulk density and heat capacity as (Global Soil Data Task Group, 2002): where ρ min [g cm −3 ] is the bulk density of the dry mineral soil. The total heat capacity of the wet soil was computed as a function of temperature and C vol (Eq. 2). Heat capacity for dry soils varied between 0.2 and 1.35 M J m −3 K −1 and averaged 0.89 M J m −3 K −1 .
The relationship between soil properties, soil water content, and the thermal conductivity λ is highly complex (Yang and Koike, 2005) and a number of empirical relationships have been developed to parameterize the thermal conductivity of different soils. A widely used relationship for the thermal conductivity of the soil as a function of soil moisture is from Hendrickx et al. (2008): where A, B, C, D are empirical coefficients that are determined based on the composition of the soil (the volume fractions of quartz, other minerals, total solids, and clay, that are, in turn, estimated from the soil fractions of sand, silt, and clay). The volume and weight fractions of sand, silt, and clay in each grid cell, are available from the standardized data set of soil horizon depth and textures for up to 15 soil horizons compiled by Webb et al. (2000). The soil depth (computed as the sum of all soil layer depths in each grid cell) in the study region generally varies between 0.0 and 9.7 m and is 4.25 m on average. Temporally and spatially varying fields of soil moisture in mineral soils were calculated using the WBM plus hydrological model (Wisser et al., 2010) that computes components of the hydrological cycle based on soil hydraulic properties and fields of climate data.
Potential evapotranspiration in WBM plus was calculated as a function of air temperature and day length using the Hamon relationship (Hamon, 1963).
We assumed that the soil is underlain by bedrock up to the lower boundary of our model (100 m), assumed the same water content as in the upper layers, and computed thermal properties like the lowest layer of the soil layers. This is similar to the approach used by  but different from the assumption of Zhang et al. (2008b) of no water in bedrock.

Thermal properties of peat soil
The heat capacity of dry peat was set to 0.58 [M J m −3 K −1 ] (Bonan, 2002), while heat capacity for wet peat was computed based on the volume fractions of peat and water and the temperature (Eq. 2). The thermal conductivities for frozen and thawed peat (λ f and λ t ) were computed using the empirical relationships found by O'Donnel et al. (2009): where θ is the water content as a fraction of saturation. The thermal properties of peat soil are therefore largely controlled by the water content in the peat soil and the position of the water table, both of which depend on the type of peatland (bog, fen), and microtopographical features. These hydrological conditions provide an important control on physical, chemical, and biological processes in peatlands (Weiss et al., 2006) and the position of the water table can both promote and hinder freezing in peatlands (Kingsbury and Moore, 1987) but very few measurements for water table or soil moisture in boreal peatlands exist (Weiss et al., 2006). In addition, the relationship between water table depth and soil moisture in the layers above the water table is site specific and difficult to generalize over a large domain.
We conceptually divided the peat soil column into fibric, mesic, and humic peat for which different levels of soil moisture θ were assumed, and the depth of those layers depends on the position of the water table in the peat column. Soil moisture for the fibric and mesic layers were set at 0.2 and 0.5 respectively and the humic peat layer was assumed to be fully saturated.
The position of the water table in peatlands is relatively stable over time (Frolking et al., 2009). Generally, the water table is deeper in bogs and higher in fens (Zoltai et al., 1998), which are the two dominant forms of peatlands in the Northern Hemisphere (Rydin and Jeglum, 2006;St-Hilaire et al., 2008). For example, Tarnocai (2006) estimates that bogs cover 67 % of the Canadian peatlands, fens cover 32 %, and swamps and marshes combined cover 1 %.
Based on observations of the mean water table from a variety of boreal peatlands from Zoltai et al. (1998) we simulated soil temperatures in peatlands for two positions of the water table; a wet scenario with a water level at 0.07 m, representing a fen, and a dry scenario with a water level at 0.30 m, representing a bog. We assumed the water table to be constant throughout the simulation period. For both scenarios, the depth of the fibric peat layer was assumed to be at 30 % of the water table depth. The different assumptions regarding the position of the water table and the depth of the soil layers are conceptually depicted in Fig. 1.

Peatland area and peat depth
To partition the grid cells into a fraction occupied by peatlands and a fraction with mineral soils, we used the recently published Northern Circumpolar Soil Carbon Database (NC-SCD; Tarnocai et al., 2009) which lists soil types and soil carbon content for some 140 000 polygons in the circumpolar region. To estimate the extent of peatland in each 30 min grid cell, we aggregated the areas of histosol and histel (frozen histosol) soil types. The resulting total peatland area for Canada and Alaska was 1.38 million km 2 in 2430 grid cells -11 % of the study area (Fig. 2). This area represents about one third to almost one half of the global peatland area north of 40 • latitude that was previously estimated to be around 3-4 million km 2 (Matthews and Fung, 1987;Aselmann and Crutzen, 1989;Yu et al., 2010).
The total soil organic carbon (SOC) in histels and histosols in the NCSCD database was 165 Pg C for the study region. The depth of these peat soils can be computed by assuming a carbon density and a peat dry bulk density, both of which vary significantly spatially depending on local conditions. For example, Turunen et al. (2001) found carbon density of 50-55 % and mean values of bulk density of 91 kg m −3 for west Siberian peatlands. Gorham (1991), based on extensive Canadian datasets, reported an average bulk density of 112 kg m −3 and a carbon density of 51.7 % of dry mass. Turunen et al. (2002) found the average dry bulk density of Finnish geological mires to be 91 kg m −3 and Vitt et al. (2009) found that peat bulk density in Canada typically varies between 90 kg m −3 and 120 kg m −3 whereas Chason and Siegel (1986) found peat bulk density to vary between 60 and 140 kg m −3 in Northern Minnesota.
To estimate the peat depth in the peat fraction of each grid cell we assumed a uniform carbon density of 50 % and peat bulk density of 100 kg m −3 . The resulting average peat depth is 3.21 m. Calculated peat depths were quite different for Alaska (generally around 0.6 m) and Canada, which had much larger spatial variability and ranged from less than 0.5 m on the Western coast to several meters in the Hudson Bay Lowlands. The total volume of peat, calculated as the product of peat depth and peat area is 3480 km 3 , 98 % of which was in Canada (Fig. 2). In cases where our estimated peat depth was smaller than 30 or 40 cm (the minimum organic layer depth for organic soils to be considered peatlands in the US and Canada), we treated those as if there was only mineral soil with a thin organic layer. This affected only 29 grid cells in Canada and Alaska (∼1 % of total number of grid cells with peatland).

Climate data
Daily time series of precipitation and air temperature for the period 2001-2100 were taken from the ECHAM 5 model (Roeckner et al., 2003) for the IPCC scenario A1B.
Simulated 20th century (20C3M) ECHAM5 air temperatures for the northern regions (20-90 • N) have been found to have the smallest bias out of all the models used in the IPCC AR4 assessment when compared to the 40 year ECMWF reanalysis (ERA40) data (Walsh et al., 2008) and was therefore chosen for the offline simulations of soil temperatures and permafrost. The original coarser resolution data was interpolated to 0.5 • grid cells using the NCAR NCL software (NCAR, 2011).
Air temperature in the northern regions is projected to increase more rapidly than in the other regions; the mean annual air temperature in the model domain increases from −3.96 • C for the first decade of the century to −1.94 • C for the period 2041-2050 and +1.0 • C for the last decade of the century. Much of the warming in mean air temperature is seen in the winter months whereas summer temperature over the entire model domain increases only slightly. Mean annual precipitation in the model domain increases from 784 mm (2001/2010) to 922 mm (2090/2100), representing an increase of nearly 18 %. This is a consistent pattern across many GCMs (Arctic Climate Impact Assessment, 2004;Meehl et al., 2007;Rawlins et al., 2010).

Model evaluation
The GIPL model has been validated against ground temperature measurements in shallow boreholes in Alaska that represent a wide range of soil and vegetation characteristics (Romanovsky and Osterkamp, 2000;Nicolsky et al., 2009). In a complementary paper (Treat et al., 2011, in prep.) the model presented here has been tested against observed time series of daily measured soil temperature at different depths in a permafrost peatland site in Northwest Territories, Canada, and a permafrost site with a thick organic horizon in Alaska (Treat et al., 2010). The results showed generally a good agreement of observed and modeled soil temperatures.
The range in seasonal mean soil temperature is generally higher in mineral soils, in particular near the surface. The differences in modeled soil temperatures are generally smaller for deeper soil layers where the temperatures are governed by long-term variations in mean annual air temperatures rather than responses to seasonal variations in air temperature.
To further evaluate the model and the assumptions regarding the properties of mineral and peat soils we compared modeled active layer thickness (ALT) under contemporary (2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010) conditions with ALT measurements from 56 sites in Alaska and Canada, representing a wide range of soil, topographical, and climatic conditions, from the Circumpolar Active Layer Monitoring (CALM; Brown et al., 2000) program. The annual ALT measurements generally represent an average value from a sample grid of 121 points over an area of 1 ha or 1 km 2 . ALT is typically measured by inserting a steel rod into the soil and recording the thawed soil depth. The maximum depth that can be measured is therefore limited by the length of the steel rod (∼1.20 m). As the period of observations in the CALM database varies for each station, we averaged annual values for the most recent 10 years and compared those values with the mean modeled annual ALT for the period 2001-2010.
ALT at a particular location depends on a number of factors including thermal properties of soil, vegetation, snow cover, soil moisture, and there is consequently a wide range of ALT observations across a spectrum of temporal and spatial scales (Brown et al., 2000) and within a single site (Hinkel and Nelson, 2003). For example, Wright et al. (2009) found the variability of the ALT in peatlands in northwestern Canada over distances of less than one meter to be correlated with lateral water flow converging to frost table depressions, a process that is not taken into account in our modeling approach. Smith at al. (2009) showed variations around the mean for one site of around 50 %. Nelson et al. (1997) reported extreme local variations in ALT, based on terrain and soil moisture conditions; for example they found a difference in ALT of 50 % between an acidic and non-acidic tundra site. Similarly, differences in the terrain lead to differences in observed ALT of 20 cm from one observation grid (1 km 2 ) to the next. ALT also varies from one year to the next.
Despite these inconsistencies in the spatial and temporal resolution of observed and modeled data that make a formal comparison difficult, our model generally reproduces the observed range of ALT well; observed average ALT ranges from 0.30 to 1.34 m (mean of 0.61 m) whereas the modeled ALT thickness ranges from 0.36 to 1.04 (mean of 0.65 m) for mineral soils, and from 0.16 to 0.42 m (mean of 0.31 m) for dry peat and from 0.18 to 0.54 (mean of 0.35 m) for wet peat. Data on soil temperatures or the composition of soil throughout the soil column are generally not available for CALM sites, so the depth of the organic layer is not known. Nonetheless, Brown et al. (2000) report that CALM data indicate that the insulating properties of organic soil generally lead to smaller interannual variations in ALT and to a shallower active layer depth, consistent with the simulations.

Trends in the timing and distribution of snow cover
Despite the projected increase in precipitation (see Sect. 2.4), the annual amount of snowfall slightly decreases over the entire model domain as a result of warming air temperatures. The mean annual values of snow depth show a significantly decreasing trend for most of the model domain: increasing trends are only seen in small regions around Northeast Alaska and in the polar regions of Canada. In areas with peatlands, the mean annual snow depth decreases on average by 0.72 mm a −1 over the simulation period.
We defined the snow season as the period when the ground is covered with a snow cover of an average depth of more than 1.0 cm over a period of five days.
The domain-mean onset of the snow season is 5 October for the period 2001-2010, but by the end of the century the domain-mean onset is 11 October. The domain-mean date for snowpack disappearance is 18 May for 2001-2010; it occurs 10 days earlier by the end of the century (8 May). Combined, the period when the ground is covered with snow decreases by around 16 days by the end of the century, representing a decrease in the snow season length of about 7 %.

Trends in soil moisture
Both the heat capacity and the thermal conductivity of soil increase with increasing soil moisture (Eqs. 5 and 6). Increasing soil moisture in the future could therefore lead to faster responses of soil temperatures to changes in air temperature (as a result of higher thermal conductivity), resulting in an increased warming of the soil. Increased heat capacity, on the other hand, leads to a slower response of soil temperatures, in particular if temperatures are low and around the freezing point (Eq. 1) and large quantities of heat need to be added/removed for the thaw/freeze phase change (Boike et al., 2009)

Trends in the frost-free season
We defined the onset of the thaw season as the first day when the soil temperature (averaged over a period of five days) at a depth of 5.0 cm rises above 0 • C and the end of the thaw season as the day when the averaged soil surface temperature drops below 0 • C. The spatial distribution of the onset of the frost-free season in the model domain follows the temperature gradient and varies from March in the southernmost parts of Canada to the end of May/early June in northern Canada. Averaged over the model domain, the frost-free period under contemporary conditions starts on 14 April and lasts until 26 September. By the end of the century, the onset will occur 25 days earlier (22 March), with the end of the frost-free period shifting from 26 September to 18 October (Table 2 and Fig. 3). The changes in the frost free season in areas covered with peatlands show a very similar pattern (Fig. 4). The differences in peat soils and mineral soils are caused both by the different thermal properties and the spatial distribution of peat soil in relation to the entire model domain (Fig. 2); peat soils are concentrated in the southern (and warmer) regions of the model domain.
The frost-free season increases from 174 to 214 days for the entire model domain, and for peatland areas from 179 to 218 days (Table 2). These increases represent a lengthening of the biologically active period of 25 % and 22 %. On average, the increase is caused by both an earlier thawing date as well as a later freeze date; the freeze date is generally more variable across the model domain (Fig. 4). The increases in the thawing period length are largest in southern Alaska, and eastern and western Canada and smallest in northern Canada, where the thawing season is shortest (Fig. 3).

Contemporary conditions
Our modeled estimate of the peatland area that is underlain by permafrost is based on the assumption of negative temperatures for at least 24 consecutive months; a grid cell can therefore be either classified as (continuous) permafrost or as non-permafrost. The estimated permafrost area varies greatly with the reference depth at which modeled temperatures are evaluated (Table 4). For a water level in peatlands of 0.30 m (the dry scenario) the area of peatlands that is underlain by permafrost varies from 43 % of the total area at 0.5 m to 72 % at 5.0 m. For the wet peat scenario (water level at 0.07 m), these estimates vary between 33 and 72 %. The higher thermal conductivity and consequently more variable active layer dynamics lead to a wider range of permafrost areas in mineral soils. The mineral soil underlain by permafrost ranges from 17 % at the surface (≤0.5 m depth) to more than 63 % at 5.0 m depth.
The spatial distribution of ALT for dry peat, wet peat and mineral soil in the model domain generally follows the climate gradient but shows some regional variation that reflects the depth of peat, the distribution of the thin organic layer on mineral soils, and the composition of the soil that impact thermal properties of soil (Fig. 5).

Future projections of the permafrost area
Warming air temperatures generally increase the active layer thickness and consequently decrease the area of mineral and peat soil underlain by permafrost depending on the reference depth used to define permafrost areas. The modeled thickness of the active layer significantly increases but is strongly dependent on the composition of the soil (peat versus mineral soil) and the moisture status in peat soils (high water  table versus low water table). Under the dry peat scenario (water level at 0.30 m), the ALT for near surface soil layers (depth < 2.0 m) increases from 0.36 m for the first decade of the century to 0.44 m at the end of the century, representing an increase of 22 %. For comparison, the mineral soil ALT in grid cell with peat increases from 0.70 to 0.93 m (32 %). Figure 5 shows the areas with an ALT of less than 2.0 m and indicates an almost complete disappearance of surface permafrost areas for mineral soil by 2100, whereas permafrost in peatlands is more persistent with smaller changes in ALT and permafrost area. The projected deepening of the ALT has large implications for the area classified as permafrost depending on the reference soil depth. Compared to contemporary conditions, the area of peatlands underlain by permafrost at 0.5 m depth declines from 43 % to 33 % of the total peatland area by the middle of the next century and to ∼14 % by 2100 (Table 4 and Fig. 6). In mineral soils, the surface permafrost (≤0.5 m depth) area decreases at a much faster rate than in peat soils, and will almost completely disappear by 2100. Near surface permafrost (2.0 m) in mineral soils decreases by ∼60 % over the century. When considering deeper soils, the permafrost in peatland areas is very stable; the permafrost area at 5 m depth is only reduced by 1 % by the end of the century whereas the mineral soil area underlain by permanently frozen soil at 5 m decreases by ∼20 % by the end of the century.

Future predictions of seasonally-thawed soil volume
We compute the maximum seasonal thawed volume of peat soil, mineral soil, and surface organic horizon each year as the product of the simulated ALT in peatlands or mineral soils and the area of peat and mineral soils. If the ALT is deeper than the peat or mineral soil or surface organic layer depths (Sect. 2.3), these depths are used rather than ALT to calculate thawed volumes. As the ALT is deeper than the thin organic layer overlying mineral soil almost everywhere under contemporary climate conditions, changes in the maximum thawed organic layer volume by 2100 are negligible ( Table 5). The maximum volume of unfrozen mineral soil, however, increases from 53 % of the total soil volume under contemporary conditions to more than 80 %, representing an increase in thawed volume of more than 50 % (Fig. 6 and Table 5).
Despite the moderate decrease in the permafrost areas in peat soils (Sect. 3.5.2), changes in the maximum thawed volume of peat as a result of a deepening active layer are also quite significant. Under contemporary conditions (2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010), 1850 km 3 of peat (45 %) are unfrozen in the dry peat scenario. This unfrozen volume of peat for the dry scenario increased to 2056 km 3 by 2050 and to 2511 km 3 (64 % of the peat volume) in the last decade of the century. The thawed volume of peat under the wet scenario is ∼5 % larger. The increase in maximum thawed volume of peat by 661 km 3 Fig. 6. Areas underlain by permafrost at different depth for dry and wet peat and mineral soil relative to 2001 (top; see Table 5 for absolute values) and time series of thawed soil volume as percentage of the total volumes (bottom; see Table 6 for absolute values) represents ∼20 % of the estimated peat volume in Canada and Alaska and can potentially have large implications for the soil carbon balance in the region.

Snow
The timing and magnitude of snowfall, the resulting snow depth, and the length of the snow season (the number of days with a significant snow pack on the ground) are some of the most important factors influencing the response of soil temperature to variations in air temperature (Stieglitz et al., 2003;Zhang et al., 2008a;Jorgenson et al., 2010;Lawrence and Slater, 2010) and therefore ultimately impact biogeochemical cycling rates in the soil (Schimel et al., 2004). Lawrence and Slater (2010) estimate that 50-100 % of the soil temperature changes by the end of the last century could be caused by variations in the snow state.
Generally, the response of soil temperatures to air temperatures is dampened the deeper the snow pack and a thinner showpack leads to colder winter soil temperatures. Longterm modeling results for Canada show a decreasing snow pack for the period 1850-2100 and suggest that the rate of increase in ground soil temperatures is 1-2 • C smaller than that of air temperature and that the rate of permafrost degradation is smaller than based on the changes in air temperature alone (Zhang et al., 2008a).
Changes in the length of the snow season can lead to an increased warming or cooling of the ground depending on when the snow season is shortened. A delayed offset of the snow season (earlier melt) increases the heat flux into the soil whereas a delay of the snow season onset can lead to both warming and cooling of the ground depending on the temporal patterns of air temperature and snowfall (Lawrence and Slater, 2010). If the delay is due to warmer air temperature and a delayed transition to below-freezing air temperatures, the result is a warming effect. If the delayed onset is primarily a result of changed precipitation patterns and less snow, the soil is not as insulated from colder air temperatures and the result will be lower soil temperatures.
Our modeled decrease in the mean annual snow depth is consistent with climate model projections across the Northern Hemisphere (Arctic Climate Impact Assessment, 2004; Zhang et al., 2008b;Lawrence and Slater, 2010) and with observed trends from measurements conducted in the last decades (Dyer and Mote, 2006). A lengthening of the snow free period between 3.1 and 6.4 days for the last 30 years has been observed using satellite data for the period 1972-2000 in North America (Dye, 2002). Our results indicate that the extension of the snow free season of 16 days by the end of the century is largely caused by an earlier melting date, resulting in less insulation of the soil from colder temperatures and consequently more heat flux into the soil. The shallower snowpack leads to less insulation of the soil from cold air during the winter months, resulting in colder soil temperatures in winter that partly offset the increases in air temperature.

Soil temperatures and active layer thickness
Warmer soil temperatures resulting in a deepening of the ALT have been widely observed in the Northern Hemisphere. Permafrost temperature at 20 m depth has increased by 1 to 2 • C in northern Eurasia during the last 30 to 35 years . This observed increase is very similar to what has been observed in Alaska and Canada Smith et al., 2010) where the warming varies between locations, but is typically from 0.5 to 2 • C. In the last 30 years, an increase of temperatures in deeper layers of permanently frozen soils observed in the Russian North and Alaska has resulted in the thawing of permafrost in natural, undisturbed conditions in areas close to the southern boundary of the permafrost zone. This warming occurred predominantly between the 1970s and 1990s.
Warming at the surface increases the ALT that controls a number of hydrological processes and determines the plant root depth, the access to nutrients in the soil, and therefore impacts biogeochemical processes in the soil . Zhang et al. (2008b) reported an increase in the ALT by 14-30 % from the 1990s to the 2090s in the permafrost region in Canada, depending on the climate forcing. For the entire Northern Hemisphere, Anisimov et al. (1997) reported an increase in ALT of 20-30 %.
These changes lead to changes in the areas classified as permafrost that is often used to report changes in the soil temperature as a result of warming air temperatures. Studies on the impact of global warming on soil conditions often report the changes in area underlain by permafrost Zhang et al., 2008b). Consistent with those projections, our results suggest that ALT in peat soils increase by 22 % and, in mineral soils in the same domain, by 32 %.

Permafrost areas and thawed volumes
To assess the contemporary geography of peatlands and permafrost areas in the model domain we used the permafrost classification of Brown et al. (2001), that subdivided the Northern permafrost region on the basis of the proportion of the ground that is actually frozen to various degrees (Hegingbottom, 2002) into continuous (>90 % permafrost), discontinuous (50-90 %), sporadic (10-50 %), and isolated (0-10 %) permafrost. Assuming the actual extent of permafrost in each class is equal to the midpoint of the reported range (for example 30 % for areas classified as sporadic), a total area of 418 300 km 2 (30 %) of peatlands is underlain by permafrost under contemporary conditions (Table 3) based on this present-day reference map (Fig. 2). Two thirds of the peatland areas are not permanently frozen under contemporary conditions. These results are comparable with our modeled results that are based on below zero temperatures at a given depth for a period of at least 2 years (Sect. 3.5.1) and are consistent with the 37 % of Canadian peatlands in permafrost areas estimated by Tarnocai (2006).
As the definition of permafrost is independent of the depth at which the layer is frozen and permafrost can be up to several hundred meters deep, a comparison of "permafrost areas" from models that simulate the temperature profile to very different depths might be misleading and approaches to determine permafrost areas based on model results consequently vary. For example, in Zhang et al. (2008b), a grid cell contained permafrost in a given year if the temperature in any layer in the soil column (modeled to a depth of 120 m) was never above freezing; the projected decline of those areas in Canada was 16-20 % by the 2090s.  followed the approach of Zhang et al. (1999) and defined permafrost as areas with a temperature below 0 • C for two or more consecutive years and projected a dramatic decline of Table 3. Contemporary geography of peatlands and permafrost extent in Alaska and Canada. Percentages relative to the total peatland area in Alaska and Canada (1.40 million km 2 ). Permafrost classification from Brown et al. (2001 ∼90 % in the "near surface" (top 3.5 m) permafrost area by the end of the century. From a biogeochemical perspective, the classification of areas as permafrost based on some (arbitrary) soil depth can be misleading as it ignores the volume of (peat) soil that is actually thawed on top of the permafrost table and could potentially interact with the atmosphere and hydrosphere, as well as the change the carbon budget of the soil should previously frozen soil become permanently thawed. This volume can be compared to the total volume and the volume of previously frozen (and mostly inert) soil volumes.
Our results indicate that despite a widespread degradation of permafrost and a rapid decrease in the area classified as permafrost, the thawed volume of peat increases by only 20 % by the end of the century (Table 5). Additional thawing of the thin organic layer is negligible because the deepening of the ALT will affect mostly deeper mineral soils while the thin organic layer thaws every season (Anisimov, 2007).

Frost free season
Between 2000 and 2100, we found an increase in frost-free season length of ∼40 days in peatlands across North America, most of which is caused by an earlier thaw in the spring.
These estimates are consistent with projected increases in the growing season and with lengthening of the growing season that has already been observed in northern latitudes. A recent increase in growing season length in northern latitudes has been found by others (Euskirchen et al., 2006;Smith et al., 2004). Euskirchen et al. (2006)

Implications for carbon cycling
The projected lengthening of the frost free period can have profound impacts on biogeochemical processes and the cycling of nitrogen (N) and carbon (C). For example, longer thaw seasons could change the N dynamics in arctic watersheds (Yano et al., 2010), while warmer winters increase the rates of N mineralization (Schimel et al., 2004). A longer growing season may increase the uptake for CO 2 (photosynthesis) but the effect on the net C balance depends on the hydrological conditions and the vegetation characteristics of the ecosystems considered (Strack et al., 2008). Lafleur and Humphreys (2008) observed higher CO 2 uptake at a low-arctic site in the Northwest Territories during a year with an early spring (snowmelt occurred 3 weeks earlier than other years) and warmer air and soil temperatures, as compared with the other years in the study. A synthesis of boreal and temperate eddy correlation sites showed that the gains in productivity due to a longer growing season were partly compensated by increased respiratory losses, and that the response of fluxes to longer growing season depends on the characteristics of the vegetation, with deciduous forests being more sensitive than evergreen forests (Piao et al., 2007). Richardson et al. (2010) show that respiration increases more than productivity when the growing season is longer in the fall, and vice-versa in the spring. Euskirchen et al. (2006) found that increases in growing season length were correlated with increases in net ecosystem productivity, and vegetation C. Often, labile C from deeper in the soil profile has been protected from decomposition by permafrost and is quickly mineralized with temperature increases following permafrost thaw (Zimov et al., 2006). Increasing soil temperatures, and the deepening of the active layer as a result of increasing air temperatures and changing snow dynamics will have implications for the cycling of carbon in peatlands and for the fluxes of carbon to the atmosphere and to the hydrosphere, as biogeochemical processes in peatlands are partly controlled by the freeze/thaw state of the (peat) soil.
The net effect of warming will depend on the balance of increased respiration and productivity Dorrepaal et al., 2009) as a result of warmer soils, and increased export of dissolved organic carbon (DOC) in rivers and streams in catchments with thawed peatlands (Frey and Smith, 2005).
Assuming that peat soils contain 50 % of carbon, a bulk density of 100 kg m −3 (Sect. 2.3), and a uniform distribution of this carbon with depth, our estimate of an additional volume of peat of 670 km 3 that would be thawed by 2100 implies that an additional 33 Pg of carbon would be thawed. Carbon concentrations in mineral soils show a much higher range than those in peat soils. The additional amount of thawed carbon from mineral soils can therefore only be estimated with large uncertainties. Turunen and Moore (2003) found mean concentrations of C in mineral subsoils beneath peat in central Finland. Slightly higher values were reported by Hossain et al. (2007) for Northern Canada.
Assuming a mean concentration of C in mineral soils of 5 kg m −3 (a mean value for mineral subsoils beneath peat in central Finland found by Turunen and Moore, 2003) , the additional volume of thawed mineral soils (10 700 km 3 ; Table 5) contains 53 Pg of carbon. An additional ∼1 Pg of C are contained in thawing mineral soils under (wet) peatland soils. The total amount of thawed carbon (87 Pg) represents ∼8 % of the entire carbon pool in the upper three meters of the soil the Northern Hemisphere, estimated to be around 1024 Pg (Tarnocai et al., 2009). If we assume the same distribution of carbon with depth and a similar response of soil temperature (∼20 % increase in the thawed peat volume) to warming in the Eurasian peatlands, an additional ∼32 Pg of the estimated 163 Pg stored in Eurasian peatlands (Tarnocai et al., 2009) could become biogeochemically active.
C processes are not only controlled by the transition of freeze/thaw but also by the temperatures alone, both in the subzero and above zero range. The warming of the soil even at temperatures below zero can have profound effects on the decomposition and production processes and thus the overall carbon balance of the peatlands (Gedney et al., 2004;Schuur et al., 2008;Carrasco et al., 2006), although the sensitivity of soil organic carbon (SOC) decomposition to soil temperature has recently been found to be scale dependent and lower on the global scale than expected from field experiments (Ise et al., 2008).
Winters and frozen-season processes are important for both the soil temperatures and the soil C balance. In addition to changes in permafrost degradation due to increased snowpack (Christensen et al., 2004;Payette et al., 2004;Stieglitz et al., 2003), Osterkamp (2005) found larger increases in permafrost soil temperatures during the winter than the summer. The duration of the thawed season at depth has significant implications for the annual C flux; Schimel et al. (2006) found that deeper, unfrozen mineral soil accounted for 50 % of the daily CO 2 production at sites in Alaska, and nearly 100 % of CO 2 production during the winter months in Alaska. Relatively large amounts of CO 2 were still produced into December, when surface soil temperatures dropped below −5 • C. The radiative forcing of peatlands can also be impacted by changes in the hydrological conditions: with wetter conditions, CH 4 emissions increase (Turetsky et al., 2002(Turetsky et al., , 2007 and are very sensitive to changes in soil temperature .

Uncertainties
Our modeling approach, like other large scale efforts to simulate soil temperatures under changing conditions (Lawrence and Slater, 2005;Zhang et al., 2008b;Wania et al., 2009) is limited by a number of simplified assumptions and the omission of small scale processes that are not represented in the model. One of these shortcomings is the neglect of convective heat transport through infiltration from rainfall and snowmelt. Although heat conduction is the dominant heat transfer process in soils (Boike et al., 2009), infiltration accelerates the warming of the soil and could become more important in a wetter climate. Infiltration is a key process in arctic soils and is controlled mostly by the thaw status of the soil , and can potentially lead to a delay in modeled spring thawing (Wania et al., 2009). The vertical transport of heat by infiltration of warmer rainwater or snowmelt is not currently considered in the model. However, experimental data suggests that the impact of soil moisture on the thermal regime is primarily through its influence on thermal conductivity and that heat advection plays a minor role (Wright et al., 2009).
The soil thermal regime is further influenced by the complex interaction between climate, topography, hydrology and vegetation, all of which cannot be adequately represented with the approach used here. For example, Jorgenson et al. (2010) reported differences in observed near-surface temperatures between boreal landscapes within similar climatic regions of ∼12 • C based on those factors; our model takes into account hydrology in a simple way but ignores slope, aspect, and vegetation.
Further biases could arise from discounting the effects of vegetation that (1) limit the amount of solar radiation hitting the soil during the summer, and (2) impact the distribution and persistence of the snow cover (Smith and Riseborough, 2002). Anticipated changes in the vegetation as a result of warming air temperatures towards shrubbier vegetation could potentially provide more shading, change snowpack and alter heat flux into the soil as well (Sturm et al., 2005).
Uncertainties in our simulation are not only caused by uncertainties in the characterization of surface and subsurface parameters  but also by the climate forcing. Predictions of soil temperature depend on the forcing climate data and are subject to uncertainties in those data; the predicted soil temperatures can vary substantially (Anisimov and Reneva, 2009). For example, Anisimov et al. (2007) found that a 0.5-1.0 • C difference the zonal mean air temperature could translate in differences in the projected permafrost area of 10-20 %, comparable to the extent of changes in the permafrost area in the next decades, although increases in temperatures, and the large scale patterns of future snow depth and precipitation are consistent across climate projections.
The estimates on peat volume are based on assumptions regarding the carbon content in peat and the bulk density of peat (see Sect. 2.3). Changing the peat bulk density from 100 to 80 increases the peat volume to 4200 km 3 whereas a peat bulk density of 120 kg m −3 would yield a volume of peat of only 2800 km 3 .

Conclusions
Consistent with observations for the recent decades and with other model simulations of soil temperatures under future climate conditions, we find a widespread degradation of permafrost in Northern regions by the end of the century. Explicitly considering the distribution and depth of Northern peatlands showed that the insulating properties of peat lead to a considerably higher persistence of permafrost in peat soil compared to mineral soils and consequently delayed degradation of permafrost in peatland areas.
Our analyses bring up the important consideration of seasonality in predicting the impacts of climate change on C cycling and storage within permafrost ecosystems. We found a 25 day earlier initiation of the spring frost-free season, and a 13 day extension of the fall frost-free period.
Despite the slower rate of soil warming in peatland areas and a slow degradation of permafrost under peat soils, a considerable volume of peat, approximately an additional 20 % of the total volume of 3480 km 3 of peat in Alaska and Canada, will be thawed by the end of the century. The potential release of carbon and the net effect of this thawing will depend on the balance between increased productivity and respiration, and will be mitigated by peat moisture.
Permafrost thaw often leads to collapse areas overlying ice-rich permafrost, and the formation of fens within existing permafrost peatlands (e.g. Beilman et al., 2001;Turetsky et al., 2007) or in areas that were formerly forested (Jorgenson et al., 2001). Jorgensen et al. (2001) estimate a 9 % increase in peatland area in the Tanana Flats of interior Alaska from formerly forested areas. Often, these thaw collapse fens are more wetter, more productive and sequester more carbon following permafrost thaw (Camill et al., 2001).
Future improvements in predictive models that help understand the role of peatlands in the global climate system and their role as a source or sink of carbon under changing climate conditions will therefore depend on coupling the thermodynamical model with a hydrological model with an appropriate parameterization of peat hydraulic properties that is able to simulate the soil moisture and water table dynamics in peatlands. As this parameterization depends on the type peatland (e.g., fen vs. bog), improved modeling approaches will require a consistent methodology to map peatlands at the global scale.