The effects of fire on the thermal stability of permafrost in lowland and upland black spruce forests of interior Alaska in a changing climate

Fire is an important factor controlling the composition and thickness of the organic layer in the black spruce forest ecosystems of interior Alaska. Fire that burns the organic layer can trigger dramatic changes in the underlying permafrost, leading to accelerated ground thawing within a relatively short time. In this study, we addressed the following questions. (1) Which factors determine post-fire ground temperature dynamics in lowland and upland black spruce forests? (2) What levels of burn severity will cause irreversible permafrost degradation in these ecosystems? We evaluated these questions in a transient modeling–sensitivity analysis framework to assess the sensitivity of permafrost to climate, burn severity, soil organic layer thickness, and soil moisture content in lowland (with thick organic layers, ∼80 cm) and upland (with thin organic layers, ∼30 cm) black spruce ecosystems. The results indicate that climate warming accompanied by fire disturbance could significantly accelerate permafrost degradation. In upland black spruce forest, permafrost could completely degrade in an 18 m soil column within 120 years of a severe fire in an unchanging climate. In contrast, in a lowland black spruce forest, permafrost is more resilient to disturbance and can persist under a combination of moderate burn severity and climate warming.


Introduction
The largest reservoir of global terrestrial carbon is stored in soils of the boreal forest (Apps et al 1993, McGuire et al Content from this work may be used under the terms of the Creative Commons Attribution 3.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI. 1995, Zoltai andMartikainen 1996, Alexeyev andBirdsey 1998). Black spruce is the dominant tree species in interior Alaska boreal forest, covering approximately 44% of the landscape (Cleve and Viereck 1983). Most of the black spruce forest in Alaska widely overlaps the discontinuous permafrost zone . Permafrost preserves the carbon in the frozen state and as a result protects it from decomposition. A thick surface layer of organic soil horizons insulates and protects permafrost from thaw (Alexeyev and  Birdsey 1982, Bonan andShugart 1989, Yoshikawa et al 2003). Therefore, the existence and thickness of the surface organic soil horizons are important factors controlling soil temperature and permafrost stability in the discontinuous permafrost region.
Climate is a major factor that directly influences the thermal stability of permafrost (Camill 2005, Callaghan et al 2011. Statistical analysis of climatic records  for Fairbanks indicate that mean annual air temperatures (MAAT) have increased by +1.79 • C during the last 80 years (figure 1(A)), and that there is no trend in snow depth over the 1930-2010 time period (figure 1(B)). Global Climate Models (GCMs) predict large temperature increases as well as increases in snow fall in high latitude regions of the Northern Hemisphere during the 21st century (Solomon et al 2007). For the moderate A1B carbon emission scenario, GCMs project annual mean temperature changes in northern high latitudes of +2.5 to +7 • C (Overland et al 2011).
In addition to the direct effect of climate on permafrost, wildfire influences the thermal state of frozen ground through burning the soil organic layer (Zhuang et al 2002, Yoshikawa et al 2003. The effects of wildfire on permafrost thermal dynamics depend on the thickness of organic soil layers remaining after fire (Yoshikawa et al 2003, Kasischke and, as a thick residual organic layer may protect underlying permafrost from thaw. The greatest risk of severe organic soil combustion occurs late in the growing season when active layer is the deepest and upper soils are dry (Stocks et al 2002. Alteration of soil organic layers in response to a changing fire regime may trigger long-term changes in permafrost stability. Although, wildfire frequency is expected to increase in Alaska during the remainder of the 21st century (Balshi et al 2009), it is possible for surface soil organic layers to recover after fire and provide necessary resilience against permafrost degradation in a warming climate.
Wildfires can alter soil moisture balance. Soil moisture content tends to increase after fire because of decreased evapotranspiration (Klock and Helvey 1976, Tiedemann et al 1979, Moore and Keeley 2000. Kasischke et al (2007) used remote sensing to observe the increase of soil moisture content in flat uplands after wildfires. In situ measurements of soil moisture made at burned sites indicate higher moisture content in comparison to unburned sites (Yoshikawa et al 2003). The increase in soil moisture content is more pronounced immediately after fire and slowly decreases within about one decade after a fire (Yoshikawa et al 2003).
Lowland black spruce forests may provide a favorable environment for the persistence of permafrost in the face of changing climate and fire regime (Robinson and Moore 2000). The low rates of evapotranspiration of lowland black spruce forest (Bonan 1991, Liu et al 2005 can cause the forest-floor to hold more moisture which favors deep accumulation of soil organic layers over time (Fenton et al 2005). The accumulation of soil organic layers in mesic-to-moist boreal forest ecosystems is associated with feedbacks among cool and moist soils, low rates of decomposition and nutrient cycling, and high moss productivity (Johnstone et al 2010a). In contrast, boreal forest ecosystems with shallow organic layers and drier soil conditions have higher rates of decomposition, higher vascular plant productivity, and lower productivity of mosses (Johnstone et al 2010a).
Our primary goal in this study is to better understand the effects and interaction among climate warming, burn severity, the post-fire recovery of surface organic soil layers, and post-fire patterns of soil moisture on the thermal stability of permafrost in lowland and upland black spruce ecosystems in interior Alaska. To assess the effect of each factor on permafrost thermal stability in black spruce forests, we conducted a series of simulations at upland and lowland sites located at the Bonanza Creek Experimental Forest located near Fairbanks, Alaska, where a forest fire occurred in the summer of 1983. In this study, we addressed the following questions: (1) which factors determine the post-fire ground temperature dynamics in lowland and upland black spruce forests?; and (2) what levels of burn severity will cause irreversible permafrost degradation in these ecosystems?

Methods
To evaluate the effect of climate warming and fire disturbance on permafrost in black spruce forests of interior Alaska we examined the following effects: (1) the effects of climate warming with no fire disturbance; (2) the effect of different levels of burn severity under a baseline unchanging climate scenario; (3) the effect of the recovery of surface soil organic layers after fire under the baseline climate scenario; (4) the effect of post-fire soil moisture dynamics under the baseline climate scenario; and (5) the effect of climate warming with fire disturbance in combination with dynamic organic soil layers and soil moisture. To simulate the effect of each these factors on permafrost we used Geophysical Institute Table 1. Thermal properties used in ground temperature simulations for the unburned upland and lowland sites. (Note: VWC-volumetric water content (fraction of 1). UWC-unfrozen water coefficients. C t /C f -thawed/frozen volumetric heat capacities (J m −3 K −1 ). k t /k f -thawed/frozen thermal conductivities (W m −1 K −1 ).) The GIPL model uses the effect of snow layer and subsurface soil thermal properties to simulate ground temperatures and active layer thickness (ALT) by solving the 1D heat diffusion equation with phase change. The phase change associated with freezing and thawing processes occurs within a range of temperatures below 0 • C, and is represented by the unfrozen water curve (Romanovsky and Osterkamp 2000). The model employs a finite difference numerical scheme over a specified domain. The soil column is divided into several layers, each with distinct thermo-physical properties. Note that GIPL is driven by a prescribed soil moisture. The GIPL model has been successfully validated using ground temperature measurements in shallow boreholes across Alaska Osterkamp 2000, Nicolsky et al 2009).
The main reason for employing the GIPL model in this analysis is that the model represents the effects of unfrozen water, which tend to be important in the modeling of the phase change processes in the frozen ground (Romanovsky and Osterkamp 2000). The GIPL model is sensitive to the changes in thermal properties of the soil layers, in particular thermal conductivity and heat capacity, both of which are depend on the unfrozen water function. The unfrozen water coefficients (UWC) shown in table 1 are dimensionless parameters obtained from calibrating the unfrozen water curve based on measurements of soil moisture. A more detailed description of the unfrozen water function calculation can be found in Romanovsky and Osterkamp (2000), Nicolsky et al (2007), and Nicolsky et al (2009). The terrestrial ecosystem model (TEM) is a processbased ecosystem model that was designed to simulate carbon and nitrogen pools of the vegetation and the soil, and carbon and nitrogen fluxes among vegetation, soil, and the atmosphere (Raich et al 1991, McGuire et al 1992. TEM represents how bio-geochemical dynamics of northern high latitude ecosystems are affected at seasonal to century scales by processes like soil thermal dynamics (Zhuang et al 2001(Zhuang et al , 2002(Zhuang et al , 2003, snow cover (Euskirchen et al 2006(Euskirchen et al , 2007, and fire (Balshi et al 2007, Yuan et al 2012. The version of TEM used in this study, DOSTEM, computes changes in organic layer thickness from calculated changes in carbon pools following fire and during stand succession. The simulation of hydrological and soil thermal dynamics by DOSTEM have been validated in Yi et al (2009b). The DOSTEM approach to simulating soil thermal dynamics differs from that of GIPL in that it calculates the heat propagation in the ground and then calculates the thaw depth (Woo et al 2004). However, DOSTEM does not consider the effects of unfrozen water on soil thermal dynamics. Soil organic thickness is computed by DOSTEM after soil carbon pools are altered by ecological processes (i.e. litterfall, decomposition) and fire disturbance based on the relationships between soil carbon content and soil organic thickness of different organic horizons in black spruce stands in Manitoba, Canada (Yi et al 2009a). This study is based on the effects of the Rosie Creek fire, a human-caused wildfire that burned 8600 acres including about one third of the Bonanza Creek Experimental Forest (Juday 2010). The Bonanza Creek Experimental Forest is located about 20 km south-west of Fairbanks, Alaska. To distinguish the effect of fire on shallow versus deep organic layers we considered upland and lowland sites being studied by the Bonanza Creek LTER Program (Jorgenson et al 2010). The organic soil profile at each site includes three organic layers (moss, dead moss, and peat) which we classify as moss, fibrous, and amorphous soil layers (table 1). The rest of the soil column includes mineral soil and rock layers. The overall soil column thickness was 20 m, where a zero heat flux was set at the lower boundary.
The baseline climate for driving the surface boundary conditions of our simulations is based on observations from 1983; the MAAT in 1983 was slightly higher than the trend line in figure 1(A) and the season snow depth was slightly below the trend line in figure 1(B). To create the baseline climate scenario, we replicated averaged daily air temperatures and snow depth from 1983 for 120 years.
The 120 year length of our simulations was chosen as a representative fire free interval for interior Alaska black spruce forest (Yarie 1981, Dyrness et al 1986, Johnstone et al 2010b. Prior to each simulation, we used pre-burn soil properties (table 1) to equilibrate initial ground temperatures by spinning-up the model until mean annual ground temperature at every depth was stabilized. The following paragraphs describe the five cases for which we conducted simulations.
(1) At the first stage of the numerical experiment we tested the effect of warming climate on permafrost without fire disturbances. In our experiment we considered a uniform linear warming trend of the MAAT by adding a yearly positive increment to the daily air temperatures within each year so that mean annual temperature increased by either 1 or 2 • C over 120 years. The air temperature warming trend was imposed gradually over the 120 years.
To introduce the effect of increase in snow fall on the thermal stability of permafrost we linearly increased snow depth by 20% over the fire frequency interval. It is important to recognize that we did not consider changes in seasonality in our simulations, such as the greater winter warming than other seasons that is generally found in scenarios based on earth system model simulations (Joyce et al 2011).
(2) To quantify the effects of different burn severities and to identify corresponding thresholds after which permafrost degradation was irreversible we removed different amounts of the surface organic layer. To identify burn severities that result in transitions from stable to unstable states of permafrost, we choose 30% (9 cm), 40% (12 cm), and 50% (15 cm) of the pre-fire organic layer thickness removals at the upland site and 25% (20 cm), 30% (24 cm), and 60% (48 cm) at the lowland site. As was stated above the most severe fire burns the amount of organic layer equal to the ALT. For example, if the overall organic layer thickness is 80 cm at the lowland site and the ALT before the fire is 48 cm, then the maximum burn severity corresponds to 48 cm of organic layer removal. Once the burn severity thresholds were identified, the next control factor was analyzed.
(3) To address the effect of new organic layer accumulation we represent the dynamics of three upper organic soil layers (moss, fibrous, and amorphous layers; table 1) as simulated by the DOSTEM (Yi et al 2009b(Yi et al , 2010. We used the DOSTEM to simulate post-fire re-accumulation rates of soil organic layers (figure 2) at the upland and lowland sites for the baseline climate scenario. The DOSTEM calculates moss layer thickness (figure 2(A)) according to an empirical function (Yi et al 2009b): where d moss is the thickness of moss (cm), d moss,max is the maximum thickness of moss (cm), y sf is number of years since fire (year), and y half = 5 is the number of years which were needed for moss to reach half of d moss,max (Yi et al 2009b). For the upland and the lowland sites we assigned d moss,max equal to 4 and 6 cm correspondingly, according to our observations. The DOSTEM simulates thicknesses of fibrous and amorphous layers (figures 2(C) and (B)) using soil carbon content, C, of each soil layer: where C is carbon content (gC cm −2 ) of an organic layer, d is organic layer thickness (cm), and a and b are fitted coefficients for the fibrous or amorphous layers (Yi et al 2009b). We integrated the organic soil thicknesses simulated by TEM ( figure 2) into the GIPL numerical model to investigate how soil organic layer re-accumulation affects the permafrost post-fire dynamics. (4) During the fourth stage of the numerical experiment we tested the effects of the post-fire soil moisture pattern on the dynamics of the active layer. The effect of the post-fire soil moisture pattern was tested together with the re-accumulation of organic layer thickness by driving GIPL with a soil moisture scenario in which we increased the upper soil layer moisture content to full saturation the year following fire and then linearly decreased the soil moisture saturation within 10 years to the pre-fire condition. This stage of the numerical experiment should be considered a first-order attempt to understand the potential sensitivity of permafrost stability to the pattern of post-fire soil moisture. (5) Finally, we tested the effect of warming climate on post-fire permafrost thermal dynamics combined with changes in soil moisture within 10 years after fire and the re-accumulation of soil organic layer thickness. Climate warming scenarios applied to the upland and lowland sites after fire were similar to the ones used at first stage of the analysis with no fire disturbances.
We conducted an independent model validation for burned and unburned black spruce sites in the general study area of our simulation sites. To our knowledge, there are no continuous measurements of the soil organic layer thickness measured at the specific site where measurements of the ground temperatures or active layer depth are available. Therefore the results of the current GIPL model validation do not include the changes in the organic layer thickness. We validate the model for the burn and unburned sites within the Bonanza Creek LTER area. Both permafrost observation stations were installed during different time periods and continuous measurements of the ground temperatures are available for the burned site from June 2000 to June 2004 and for the unburned site from March 1996 to August 1999. These permafrost observation stations represent a small climate station which includes high-precision air and ground temperature sensors (Campbell Scientific L107 thermistor and MRC multi-thermistor probe) and up to three Hydra Probe soil moisture sensors. Ground temperatures are measured down to the 1 m depth with sensors located about every 0.1 m. All measurements are taken at 1 h time intervals. The thermistor sensors were calibrated in an ice bath prior to installation to an accuracy of 0.01 • C. More detailed information on data logger installation can be found in Osterkamp and Romanovsky (1999), Romanovsky and Osterkamp (2001), Osterkamp (2003), and Romanovsky et al (2003).
The burned validation site is located within the Bonanza Creek LTER site on the Tanana River floodplain, 26 km south-west from Fairbanks. This site experienced a severe forest fire in 1983 during which the entire organic layer was burned. The calculated thaw depth during the 2000-2004 time period was about 4 m thick which correlates well with our observations. Current vegetation is young small black spruce trees with shrubs. The results of burned site ground temperature simulations in comparison with measured data are shown in figure 3. To quantify the difference between measured and simulated ground temperatures we calculated the mean average error (MAE) at every depth where temperatures were measured. The MAE represents a sum of the absolute daily differences between measured and simulated ground temperatures over all of the days contributing to the average.
The unburned validation site is located within a black spruce stand on the old floodplain of the Tanana River about 1.5 km south-west from the burned site. Vegetation includes sparse black spruce stands with shrubs and a thick moss cover dominated by Sphagnum. The calculated active layer thickness during the 1996-1999 was about 0.6 cm thick. The results of unburned site ground temperature simulations are shown in figure 4.

Results
The results of ALT simulations under different climate warming scenarios without fire disturbance or changes in snow thickness indicate gradual thickening of the active layer for every warming scenario at the upland and lowland black spruce sites (figure 5). For the upland simulation with increased snow thickness (+2 • C MAAT warming over 120 years and 20% increase in snow thickness) there was rapid increase in permafrost degradation at approximately year 105 when snow increased by 17.5% (see figure 5(A)). The thickening of the active layer was weaker at the lowland site for all warming scenarios, and there was no rapid increase in permafrost degradation for the scenario with increasing snow thickness ( figure 5(B)).
Fire disturbance with no climate change and no organic layer regrowth had a substantial impact on permafrost thermal stability. In the upland black spruce site permafrost started to thaw progressively when 12 cm (30%) of the pre-fire organic layer thickness was removed (figure 6(A)). In contrast, permafrost did not progressively thaw at the lowland site until 24 cm (30%) of the pre-fire organic layer thickness was removed (figure 6(B)). Fire disturbance at the lowland site did not immediately degrade permafrost after 24 cm of the pre-fire organic layer was burned until approximately 100 years had passed. It is notable that for the highest burn severity (48 cm, 60% removal) at the lowland site the rate of permafrost degradation was less pronounced than for 15 cm (50%) removal at the upland site.
The next step in the numerical experiment was the implementation of the dynamic change of the organic soil layer horizons in the GIPL model after fire, under the baseline climate scenario. When the organic soil layers were allowed to re-accumulate after fire at the upland site, the permafrost did not totally degrade down to 5 m for 15 cm (50%) removal of pre-fire organic layer thickness ( figure 7(A)). In contrast, permafrost at the lowland site was able to fully recover its thermal state and ALT within 40 years after fire for all of the levels or organic layer removal we considered ( figure 7(B)).  In comparison to the reference simulation, the simulation with post-fire soil moisture dynamics had a subtle effect on the pattern of ALT dynamics (figure 8). A talik, the unfrozen ground layer that survives winter, formed couple years earlier at the upland permafrost site for the simulation with dynamic soil moisture, contributing to a slightly deeper permafrost table ( figure 8(A)). In contrast, at the lowland site the talik did not form and an excess of soil moisture during 10 years after fire slowed the thickening of the active layer ( figure 8(B)).  More detailed analysis of the dynamic soil moisture effect indicates that the increase in soil moisture content reduced the ground temperature seasonal variability during the year following fire ( figure 9(A)) in comparison to no change in soil moisture ( figure 9(B)). In the current formulation of GIPL, changes in the soil moisture influence permafrost dynamics through the effects of the latent heat of water fusion, but not soil thermal properties (Jorgenson et al 2010). The changes in latent heat after fire corresponds to the energy required to thaw or freeze the excess of ground water available as a result of a decrease in evapotranspiration. Thus this latent heat effect dampens the ground temperature seasonal variability ( figure 9(A)). At the upland site with 15 cm organic layer removal, changes in soil moisture content accelerate development of the thaw layer and therefore increase the warming effect of the ground. The increase in soil moisture affects the active layer depth and requires more time for the ground to refreeze. The analysis of simulated ground temperature freeze-up date for increased post-fire soil moisture content indicates a faster transition from seasonally thawed active layer to talik (figure 10(A)). Four years after fire the freeze-up is not happening anymore for the dynamic moisture simulation. This threshold marks the time of the beginning of talik formation, which allows heat to stay longer in the ground and contributes to its further warming. In the case of static soil moisture it takes 2 years longer for the talik to form (figure 10(B)) which reduces the effect of the ground warming and as a result reduces the depth of the permafrost table ( figure 8(A)). At the lowland site the increase in the soil moisture delayed the freeze-up date only slightly and no talik formed ( figure 10(B)). The ability of the lowland active layer to completely refreeze and the increased ice content in this layer contributed to the shallow ALT within 10 years after fire ( figure 8(B)).
The simulations that combined climate warming and wildfire including all the post-fire effects described above indicate substantial permafrost degradation at the upland site for 15 cm of organic layer removal (figure 11(A)) in comparison to no organic layer removal ( figure 5(A)). In contrast, permafrost at the lowland site for 48 cm of organic layer removal was able to recover its pre-fire thermal condition under almost every warming scenario, except the highest warming scenario with an increase in snow thickness ( figure 11(B)). For the scenario of 100% removal of the pre-fire organic layer at the upland site, the ALT dramatically increased to 8 m within 30 years after fire for the no

Discussion
Resilience of permafrost is the capacity to maintain temperatures below freezing and stable ground ice contents and morphologies, whereas vulnerability is the extent to which permafrost thaws vertically and laterally and how much thaw settlement occurs during thawing of ground ice (Jorgenson et al 2010). In the current analysis we defined the vulnerability and resilience of the permafrost based primarily on changes in its temperatures and ALT. Permafrost vulnerability in black spruce forest is the combination of resilience and exposure, where climate change and burn severity are components of exposure.
The magnitude of the effects of climate warming and burn severity on permafrost can be influenced by soil texture, thermal properties, the degree of soil saturation, and snow depth. An increase in snow depth enhances the effect of climate warming by insulating the ground from winter air temperatures. However, the overall effect of snow after fire is not well understood. For example, snow compaction in young forests may decrease the insulative effectiveness of snow if open regrowing forests are more subject to the effects of wind after fire. In contrast, forests intercept more falling snow as they age, which may improve permafrost stability in old forest due to less ground insulation.
The ability of permafrost to recover after fire depends on the thickness of the organic layer remaining after fire. Our simulations indicate that the effect of the re-accumulation of organic layer thickness after fire increased permafrost resilience for both boreal forest sites. Furthermore, our analysis indicates that recovery of the permafrost thermal conditions after moderate to severe fires could be initiated as soon as moss layer starts to re-accumulate. Thus, the amount of organic layer left after fire and the re-accumulation of the soil organic layer are two major negative feedbacks that provide resistance against permafrost degradation.
Fire is a major disturbance in boreal forest that influences the vegetation and permafrost (Balshi et al 2009). Wildfires are strongly influenced by climate and human activity (Kasischke et al 2000). Empirical models of fire frequency and severity in Alaska and Canada suggest that average area burned per decade will double by 2041-2050 and will increase on the order of 3.5-5.5 times by 2100 (Balshi et al 2009). The combustion of vegetation, moss, and a portion of the surface peat layer by fire led to an increase of the ALT and degradation of permafrost. For ecosystems with thin soil surface organic layers, the degradation of permafrost could be rapid and irreversible in a warming climate. Permafrost degradation could also be affected by patterns of ecosystem succession after fire depending on whether black spruce forests tend to replace themselves after fire disturbance or whether they transition to deciduous forests with thinner organic layers than black spruce forests Kasischke 2005, Johnstone et al 2010b). Furthermore, the thawing of permafrost increases the depth of active layer and lowers of the water table, which exposes carbon in permafrost soils to decomposition (Schuur et al 2009). Other studies have noted that the combination of climate warming and increased fire frequency have the potential to release carbon from boreal forest soils in amounts that might have consequences for efforts to mitigate the buildup of greenhouse gases in the atmosphere (Hayes et al 2011, Yuan et al 2012. Through studying the resilience of permafrost in response to climate warming and fire, our study further elucidates why it is important to represent interactions among climate warming, fire severity, permafrost dynamics, and soil carbon dynamics of boreal forests in earth system models that are used to assess how responses of boreal forests to climate change may influence the climate system.

Conclusions
The thawing of permafrost can have a significant impact on soil hydrology, vegetation succession and on the global carbon balance. The results of the simulations in this study indicate substantial vulnerability of discontinuous permafrost in the black spruce forest with thin (up to 30 cm) organic soil layers.
Our analysis experiments showed that soil organic layer dynamics plays a crucial role in permafrost recovery after fire disturbances. The predicted increase in forest fire frequency and severity with climate warming will most likely accelerate permafrost degradation in forests with thin and relatively dry organic soil layers. The thickness of soil organic layers left after fire and the rate of its recovery will determine the rate of permafrost degradation and whether or not permafrost will recover, partially degrade, or fully degrade. Wet soils provide a better environment for the regrowth of the moss layer after fire, and therefore promote permafrost resilience. However, the increase in the post-fire soil moisture content and changes in soil thermal properties have the potential to accelerate permafrost degradation through effects on heat transfer in the soil. In the version of the GIPL model used in this study, the effects of changes in soil moisture were not fully addressed, since changes in soil moisture were not coupled with changes in soil thermal properties. Future work is needed to better address changes in soil moisture content in the GIPL model. In this study we primarily addressed how increases in soil moisture influence permafrost dynamics through latent heat effects. This study emphasizes the importance of interactions among climate warming, soil organic layer thickness, and fire severity in influencing the vulnerability of permafrost in boreal forests. These interactions are important to consider in assessing the response of carbon dynamics in boreal forests to projected changes in climate and fire disturbance.