Long-term ecological consequences of forest fires in the continuous permafrost zone of Siberia

Wildfires are an important factor in controlling forest ecosystem dynamics across the circumpolar boreal zone. An improved understanding of their direct and indirect, short- to long-term impacts on vegetation cover and permafrost–vegetation coupling is particularly important to predict changes in carbon, nutrient and water cycles under projected climate warming. Here, we apply dendrochronological techniques on a multi-parameter dataset to reconstruct the effect of wildfires on tree growth and seasonal permafrost thaw depth in Central Siberia. Based on annually-resolved and absolutely dated information from 19 Gmelin larch (Larix gmelinii (Rupr.) Rupr.) trees and active soil layer thickness measurements, we find substantial stand-level die-off, as well as the removal of ground vegetation and the organic layer following a major wildfire in 1896. Reduced stem growth coincides with increased δ13C in the cellulose of the surviving trees during the first decade after the wildfire, when stomatal conductance was reduced. The next six to seven decades are characterized by increased permafrost active soil layer thickness. During this period of post-wildfire ecosystem recovery, enhanced tree growth together with positive δ13C and negative δ18O trends are indicative of higher rates of photosynthesis and improved water supply. Afterwards, a thinner active soil layer leads to reduced growth because tree physiological processes become limited by summer temperature and water availability. Revealing long-term effects of forest fires on active soil layer thickness, ground vegetation composition and tree growth, this study demonstrates the importance of complex vegetation–permafrost interactions that modify the trajectory of post-fire forest recovery across much of the circumpolar boreal zone. To further quantify the influence of boreal wildfires on large-scale carbon cycle dynamics, future work should consider a wide range of tree species from different habitats in the high-northern latitudes.


Introduction
. Despite an apparent decline in the total global area burned between 1996 and 2015 (Doerr andSantín 2016, Andela et al 2017), the frequency and intensity of boreal wildfires in Alaska, Canada and Russia have increased substantially over the 20th and early-21st century (Soja et al 2007, Turetsky et al 2011, Ponomarev et al 2016, Forkel et al 2019. This evidence is in line with global and regional model output that predicts a further increase in the extent and severity of boreal wildfires due to climate change (Chapin et al 2000, Flannigan et al 2013, Boulanger et al 2014. The spatiotemporal distribution of forest fires in Siberia also testifies an increasing danger under future surface warming of the highnorthern latitudes (Ponomarev et al 2016, García-Lázaro et al 2018. Up to 80% of the boreal forest grow on permafrost, where only a shallow upper soil layer thaws temporally each summer (Cable et al 2014;Helbig et al 2016), the so-called active soil layer (ASL). Although the global permafrost extent is largely controlled by surface air temperature (Shur and Jorgenson 2007), wildfires can affect the ASL by removing the insulating upper vegetation and organic soil layer, thus facilitating vertical heat transfer (Jafarov et al 2013). To understand the effects of wildfires on permafrost, the ASL thickness was repeatedly measured at many sites in northern North America (MacKay 1995, Brown et al 2000, Viereck et al 2008, and upper permafrost thaw dynamics were estimated from the rate of forest ecosystem recovery after burning (Brown et al 2015). Previous studies demonstrate not only high spatial variability, but also great dependency of the ecological consequences of wildfires on a multitude of environmental factors, including seasonal permafrost thaw, soil texture and moisture, as well as the timing and intensity of fires (Minsley et al 2016). Despite the above, our understanding of the longevity of postwildfire ecosystem recovery is still limited (Shvetsov et al 2019), because comprehensive and interdisciplinary long-term monitoring studies in the boreal forest are logistically challenging.
Dendroecology, however, can provide insights of unique temporal resolution, because tree rings may allow fire histories to be reconstructed (McBride 1983, Stivrins et al 2019. A recent example of the successful utilization of tree rings in wildfire dendroecology is the precise dating of moss buried stems to quantify postwildfire dynamics of the ASL thickness and ground vegetation recovery in northern Siberia (Knorre et al 2019). Moreover, innovate dendrochronological approaches have combined annually-resolved and absolutely-dated ring width measurements with stable isotopic ratios to provide eco-physiological insights into tree-fire interactions (Beghin et al 2012, Battipaglia et al 2014. Stable carbon and oxygen isotopes in wood cellulose can reflect information on the wateruse, stomatal conductance and photosynthesis of trees (McCarroll and Loader 2004). Despite a large body of tree-ring stable isotope research across the boreal permafrost zone (Saurer et al 2004, Kirdyanov et al 2008, Sidorova et al 2009, Tei et al 2013, Churakova (Sidorova) et al, 2019, only a few studies addressed the impact of wildfires on tree growth (Porter et al 2009, Sidorova et al 2009, and none of them assessed the direct and indirect effects of post-fire ecosystem recovery on the isotopic composition in tree rings. Here, we combine dendrochronological and stable carbon and oxygen isotopic measurements to reconstruct the impact of forest fire and climate on the radial growth and physiology of Gmelin larch (Larix gmelinii (Rupr.) Rupr.) in the continuous permafrost zone of Central Siberia. We therefore identify wildfire-driven, interannual to multi-decadal changes in various treering parameters during ecosystem recovery. We also reconstruct the prior-and post-fire dynamics in different ecosystem components, and link the observed changes in tree physiology to the recovery rate of the permafrost ASL after wildfires.

Materials and methods
Our sampling site is located within the continuous permafrost zone in the northern part of Central Siberia (64°13′N, 100°28′E and 215 m asl). This region is characterized by continental climate with a distinct intra-annual temperature amplitude between the warmest (T July =+16.6°C) and coldest (T January =−35.9°C) months, and overall very low annual precipitation totals ∼360 mm (calculated since 1929 for the meteorological station in Tura that is located about 13 km away from our sampling site). Permafrost thickness varies between 220 and 500 m (Brown et al 1997), with up to two meters of seasonally thawing ASL. The natural forest is dominated by Gmelin larch (Larix gmelinii (Rupr.) Rupr.), which is well-adapted to the harsh environmental conditions of Siberia's boreal zone (Abaimov et al 1997). The growing season is restricted to ∼70-90 d between late-May and early-September (Bryukhanova et al 2013, Shishov et al 2016. Except for a small tree island, a massive wildfire in 1896 removed most of the forest cover in our study area ( figure 1(a)). Due to an extent area of >50 km 2 that was burned, and a high rate of tree mortality, the 1896 wildfire clearly exceeded most of the recent fires that affect ∼20 km 2 , on average Ponomarev 2017, Ponomarev andPonomareva 2018). Although trees that survived the wildfire are larger than those that established afterwards ( figure 1(b)), the post-fire stand is much denser (table 1). A higher proportion of lichens, including several species of Cladonia and Cetraria genera suggests drier conditions at the older stand compared to the younger post-fire stand that is mainly covered by mosses, e.g. Pleurozium schreberi and Hylocomium splendens (Vodop'yanova 1976).
To assess the regional-specific rate of wildfireinduced changes in ASL thickness, we performed field measurements of the seasonal depth of permafrost thaw at eight sites that have been affected by wildfires in different years. While the most recent fire happened just a few weeks before the soil measurements, previous fires occurred 1, 11, 15, 25, 60, 107 and 180 years ago. All of the study sites are characterized by a specific permafrost microtopography termed 'earth hummocks' that is typical for the area. Such microrelief patterns consist of mounds neighboring bowl-shaped depressions that are formed by the freeze-thaw cycles within the ASL (Kawahigashi et al 2011 and references therein). The ASL thickness was measured in the second half of the growing season of 2006 separately for troughs and mounds along a 10 m transect within every plot at 1 m intervals (n=11 positions in total) with a steel pole (1 cm measurement interval, 1.2 m length and 1 cm diameter). To confirm our data, we additionally measured ASL thickness for randomly selected 7 mounds and 7 troughs across each study site. The ASL thickness was obtained as the difference between the measured total thawed depth and the height of the organic layer (moss-lichen stratum and organic soil layer accumulated on a mineral soil surface). The ASL thickness therefore only refers to the mineral soil. For the same sites and transect points, we measured the moisture content at 5 cm depth mineral soil. Samples were weighed, and subsamples in triplicate were dried in an oven at 105°C for 48 h to determine field moisture content as percentage of the absolute dry weight. Although these additional measurements were restricted to mid-July and mid-August, they are representative of the warm season, because the fastest rate of annual thawing of the upper permafrost layer happens before mid-July Three to four increment cores per tree were collected at the main sampling site in the summer of 2011 from eight dominant trees that survived the wildfire in 1896 (Kirdyanov et al 2008(Kirdyanov et al , 2018. Two larches were felled afterwards and cross-sections were taken at 0.3 m stem height to date the wildfire with fire scars. For the larch stand that recruited after the wildfire, we collected two cores from eleven randomly selected dominant trees. One core per tree was used for x-ray densitometry with DENDRO-2003(Schweingruber 1988. The resulting tree-ring width (TRW) and maximum latewood density (MXD) measurements were simultaneously used for visual cross-dating, which was statistically verified by COFECHA (Version 6.02 P). The individual TRW series were standardized with  negative exponential functions, whereas cubic smoothing splines with 50% frequency-response cutoff at 2/3 of the individual series length were used for the age-trend removal in MXD data. Bi-weight robust means of the individual measurement series were calculated to produce dimensionless TRW and MXD index chronologies. Moving-window, inter-series correlation coefficients (Rbar) and the Expressed Population Signal (EPS; Wigley et al 1984) were used to characterize the general behavior and internal signal strength of the individual chronologies.
To assess changes in wood accumulation rates at the two-dimensional stem surface level, the raw TRW measurements were converted into basal area increment (BAI), and the resulting BAI series of the individual trees were averaged separately for the prior-and post-wildfire cohorts. BAI was used in addition to the raw and standardized TRW and MXD chronologies, because it is less prone to age-related changes in tree geometry (Fritts 1976), and therefore does not require standardization/detrending (Biondi and Qaedan 2008).
For the stable carbon and oxygen isotope analysis, we used material from two cores of five trees that survived the fire. After cross-dating, wood of each calendar year of the individual tree-ring cores was split, and pooled in accordance with the weight contribution of each sample before milling. Further procedures of cellulose extraction, sample preparation and isotopic composition measurement followed the steps described in Saurer et al (1997Saurer et al ( , 2016. The isotopic values (δ 18 O and δ 13 C) in cellulose are expressed in delta notation relative to the international V-PDB reference for carbon and V-SMOW for oxygen. The raw carbon isotope data were corrected for the decreasing atmospheric δ 13 C value due to fossil fuel combustion (Suess Effect, Keeling 1979). Based on δ 13 C values, we calculated intrinsic water use efficiency (iWUE) according to Farquhar and Richards (1984). Further, iWUE relations to BAI were used to assess how trees are adapting to changing environmental factors after the occurrence of wildfires.
To determine the most important climatic factors that control radial tree growth and isotopic composition, the different tree-ring parameters (TRW, MXD, BAI and δ 13 C and δ 18 O) were correlated against monthly temperature means and precipitation totals from the nearest meteorological station in Tura that is operating since 1929 (WMO_ID=24507). June-August summer temperature means (Ts) and winter precipitation totals from previous November to February of a current year (Pw) were also used. To observe growth-climate response shifts, Pearson's correlation coefficients were calculated for two similarly long, non-overlapping periods : 1936-1964 and 1980-2009. While both intervals are characterized by good data quality and quantity, they represent intervals with different tree-ring growth rates. To analyze differences among mean values of iWUE and BAI for various periods of tree growth, we applied a variance analysis (one-way ANOVA). The correlation analysis, the Shapiro-Wilk W test for normality prior to ANOVA and, the ANOVA itself were carried out in Statistica, version 10 (StatSoft, Inc.).

Results
Radial tree growth Tree-ring data indicate that a devastating wildfire in 1896 caused a pulse of larch regeneration with trees reaching 1.2 m vertical height (equal to sampling height) in 1908-1927 ( figure 2(a)). Both prior-and post-wildfire trees show enhanced juvenile growth which is removed by standardization (figures 2(b) and S1, available online at stacks.iop.org/ERL/15/ 034061/mmedia). TRW index chronologies are characterized by high synchrony of year-to-year and multidecadal fluctuations, both for prior-and post-wildfire tree generations, with high inter-series correlation (Rbar=0.504, EPS>0.85 since 1873 and Rbar=0.460, EPS>0.85 since 1936) (figures 2(b), S1).
An abrupt decline in tree growth after the wildfire is well pronounced in the raw and standardized TRW chronologies. Tree radial growth indices during the pre-wildfire period (0. 95±0.24, period 1923-1895) considerably decrease to 0.54±0.12 in 1896-1907 with the minimum value of 0.31 in 1900, which is 38% compared to the last prior-wildfire year (0.82 in 1895). Between 1908 and 1973, TRW indices vary around the mean of 1.20±0.27, with a general upward trend during the first three decades of post-wildfire recovery and a continuous decline from the mid-1950s to the 1970s, when standardized tree growth drops below the mean value (0.72±0.18). Similar tendencies are observed for wildfire-induced changes in BAI (figure 2(c)). Enhanced wood accumulation rates during the last two decades before the wildfire are followed by an approximately decade-long reduction in BAI and a continuous recovery from around 1908 to the 1930s. Then, the increasing trend is switched to a decreasing tendency until the late-1970s. In the 1980s and onwards, BAI fluctuates around and below 100 mm 2 year −1 per tree ( figure 2(c)). Mean BAI of the younger post-wildfire tree generation follows the multi-decadal trends of the older trees, but is generally characterized by three-fold lower values .
MXD has been affected by the wildfire similarly to TRW and BAI ( figure 3(a)

Stable isotope ratios
The mean value of stable carbon isotope ratios in treering cellulose before the wildfire was −23.75±0.43‰ (1860-1895 period) ( figure 3(b)). The wildfire event caused an increase of δ 13 C from −24.16‰ in the pre-wildfire year 1895 to −21.87‰ in 1897, and −22.41‰ in 1898. In the following decades, there is an increasing trend of C-isotope values (0.11‰ per decade, p<0.0001) from 1897 to 1973, with the mean value of −23.11±0.57‰ overall. The 1974-2009 period starts with a rather abrupt decline of δ 13 C to the mean −23.32±0.45‰, which is close to the pre-wildfire values.
The relation of BAI to intrinsic water use efficiency (iWUE) is shown for three periods: prior-wildfire (1860-1895), period with high radial growth  and with reduced growth (1980-2009) (figure 4). The prior-and post-fire periods are similar in the rates of wood accumulation (one-way ANOVA p=0.316), but differ significantly in iWUE (p<0.01). The period from 1980 to 2009 is characterized by increasing iWUE and significantly reduced BAI (p<0.01).

Growth-climate sensitivity
Temperature effect: Correlations between the treering parameters and monthly climate data indicate a rather low dependency of tree radial growth on temperature during the 1936-1964 period of enhanced wood production with June thermal conditions affecting TRW (r=0.39, p<0.05), and previous year December (r=−0.45, p<0.05) and current January (r=0.39, p<0.05) temperature influencing δ 13 C and δ 18 O, respectively (table S1). During the 1980-2009 period of low tree radial growth, the role of temperature, especially during summer, increases. All the tree-ring parameters are positively related to mean summer temperature, with June being especially important for TRW (r=0.42, p<0.05), MXD (r=0.40, p<0.05) and BAI (r=0.45, p<0.05) and July for δ 13 C (0.52, p<0.01). In contrast, prior-season temperature in April has a negative effect on TRW, MXD and BAI (up to r=−0.51, p<0.01).
Precipitation effects: During the 1936-1964 period, TRW and BAI are negatively influenced by preseason precipitation (up to r=−0.56, p<0.01), while MXD depends positively on previous October (p<0.05) and current September (r=0.48, p<0.01) monthly precipitation. Prior-season precipitation is also important for δ 13 C and δ 18 O. The oxygen isotope ratios are negatively affected by winter (r=−0.65, p<0.01) (mostly, December and February) and July (p<0.05) precipitation while May precipitation positively influences the carbon and oxygen isotope ratios (p<0.05). July precipitation is also important for δ 18 O (p<0.05). During 1980-2009, the influence of pre-season precipitation on TRW, MXD, BAI, and δ 18 O decreases, but the positive role of July and previous October precipitation increases (p<0.05). For δ 13 C, the negative influence of July precipitation remains similar for this period but a positive influence on May precipitation and the role of accumulated winter snow amount for δ 18 O decreases.

ASL dynamics
Data on the depth of seasonal permafrost thaw for a chronosequence of forest ecosystems affected by wildfire from 0 to 180 years ago show a large and significant increase of ASL depth due to burning of the insulating layer of ground vegetation and organic layer ( figure 5(a)). Seasonal thaw depth increases from circa 30 to 50 cm to maximum 1.5-1.7 m within the first 30 years after the fire. Gradually developing ground vegetation covers the soil surface approximately 20 years after the wildfire and continue to accumulate enhancing the insulating effect. Therefore, ASL depth decreases and returns to the prior-wildfire depth, but only 70-80 years after fire. Furthermore, there is great difference in ASL thickness between mounds and troughs with the former showing greater permafrost thaw up to 157±14 cm compared to troughs with 103±18 cm. Moisture of mineral soil at 5 cm depth decreases to a minimum within the first 20-30 year after wildfire and then slowly recovers ( figure 5(b)). Mounds are generally drier and their soil moisture decreases significantly faster after the wildfire than in troughs.

Discussion and conclusions
Our evidence for post-wildfire dynamics of the ASL thickness and changes in tree-ring parameters reveal a three-phase recovery process of Siberian forest ecosystems after burning ( figure 5). During the first phase, 10-15 years following a wildfire in 1896, TRW, BAI, and MXD declined, while δ 13 C was higher than before the fire. In this period of reduced radial tree growth that is longer than previously reported for the area located slightly north (Kharuk et al 2011), the negative effects of the fire on growth and physiology appear to be dominant. Though most of Siberia's high-latitude wildfires are ground fires (Conard et al 2002, Sofronov andVolokitina 2010), there is a direct thermal damage on larch needles. In addition, even a low intensity wildfire can cause serious damage to trees on permafrost because of the superficial root system within the soil organic layer (Kajimoto et al 2003). Direct fire damage of larch roots may lead to a complete destruction of the tree stand over large territories and severe disturbance for a few surviving trees. The immediate post-fire increase of δ 13 C ( figure 3(b)) reflects lower leaf-internal CO 2 -concentration, higher water-use efficiency, and hence less tree discrimination against heavier carbon isotopes, which could be caused by stomata closure or faster tree photosynthesis rates (McCarroll and Loader 2004). However, low tree-ring growth (BAI) and MXD, partly depending on the amount of carbohydrates accumulated during the growing season (Rinne et al 2015b), testify against enhanced photosynthesis during the first years after fire. Therefore, the assumption of wildfire-caused damage to tree roots impeding the water supply supports the hypothesis of stomata closure and less discrimination against 13 C by Rubisco. The increased values of carbon isotope ratios in tree rings during the first years after the wildfire in 1896 can also be connected to the enhanced use by trees reserve carbohydrates with a heavier carbon isotopic composition due to the damage of the photosynthetic apparatus (Rinne et al 2015a).
The second phase from approximately 1908 onwards covers the following six to seven decades and coincides with the period of increased ASL thickness and elevated tree growth (figure 5). Thermal soil conditions are of great importance for tree growth on permafrost (Nikolaev et al 2011, Bryukhanova et al 2013, Churakova (Sidorova) et al, 2016, and deeper ASL usually promotes faster tree growth due to higher soil temperatures near the surface and higher water and nutrient availability , Bryukhanova et al 2015, Shishov et al 2016. At our study site, a favorable hydrothermal regime of soil stimulated successful regeneration of larch, development of deeper rooting systems of the survived trees, and higher tree radial growth. Interesting, BAI increased during the first half of the period till the 1930s and gradually decreased afterwards, which generally reflects the ASL dynamic ( figure 5(a)). Climate correlation analysis testifies that during this period, summer temperatures only slightly influence TRW and have no effect on other tree-ring parameters which implies favorable soil temperature regime for tree growth.
The positive trend in tree-ring δ 13 C to heavier values combined with a declining trend of δ 18 O during the period of deeper ASL may reflect the increasing photosynthesis rate with little changes in stomata conductance when assuming that mainly leaf-level isotope fractionations dominate the isotope signals (Scheidegger et al 2000). The increased intensity of photosynthetic activity is supported by higher wood accumulation (figure 2), which is in the range of priorfire values though intrinsic water use efficiency is significantly higher (figure 4). Increasing carbon isotopic values in the tree-rings could also be influenced by soil-respired CO 2 that is re-assimilated by larch trees (Buchmann et al 1998). Respired CO 2 is rather depleted compared to atmospheric CO 2, and therefore changing uptake over time may result in assimilated carbon trends. However, this effect was found to be significant close to the soil and should therefore not influence the isotope ratio of larch trees after the juvenile phase. Due to the fact that δ 18 O values of plant tissues partly reflect isotope composition of source water (Roden et al 2000), the interpretation of the tree-ring oxygen isotopic ratios in the permafrost zone is complicated by the ability of trees to utilize lighter water from thawing snow and permafrost (Sugimoto et al 2002, Saurer et al 2016. The declining trend of δ 18 O may also indicate an increasing role of meltwater supply for trees. This is further confirmed by a strong negative relation between δ 18 O and winter precipitation during the period of deep ASL and enhanced tree growth (table S1). An increased variability of tree-ring oxygen stable isotope ratios during the first years of post-fire recovery from 1904 to 1920 could reflect varying amounts of depleted melted versus enriched rain water (Saurer et al 2016) while undeveloped and patched ground vegetation cover still recovering after the wildfire is not able to intercept the rain water. During this period, the water supply of trees therefore mainly depends on summer conditions, when the lack of water during dry summers may be compensated by a higher use of permafrost water (table S1). Importantly, tree growth and isotopic composition depend much on winter and pre-growing season precipitation. Though soil moisture is relatively low ( figure 5(b)), the influence of summer precipitation on tree growth is not pronounced, which implies an easier access of atmospheric water to soil and tree roots.
The third phase, the period of reduced tree radial growth from the 1970s, is characterized by a thinner ASL due to the recovered ground vegetation, and a higher summer climate limitation of tree growth (table S1). Low carbon isotope values associated with low BAI (figure 4) testify the decreased photosynthesis rate that is caused by low soil temperature and obstructed water supply due to shallow ASL and low soil temperatures ( figure 5). An abrupt decline of δ 13 C in the early 1970s is probably associated with a tipping point in the adaptation of larch trees to the soil conditions that become similar to those in the prior-wildfire period. This fast change could, for instance, be related to the die-off of parts of the roots due to shallower ASL (Knorre et al 2019). Increased iWUE and reduced wood accumulation (BAI) (figure 4) demonstrate a decline in photosynthetic rate as a result of adaptation to cooler soil conditions. For this period, we assume a complete recovery of the forest ecosystem after the wildfire. BAI of younger post-fire trees is, however, still lower than basal growth of the prior-wildfire cohort, both recent and at the similar ages of priorwildfire trees in the past (figure 3), evidencing an even longer effect of wildfires on forest carbon accumulation. On the other hand, low carbon accumulation rates by single trees are compensated by higher stand density (table 1), which confirms the importance of a trajectory for the post-wildfire ecosystem recovery (Alexander et al 2012(Alexander et al , 2018. The results of our dendroclimatic analysis show temporal changes of tree growth response to climatic variables (figures 4, 5), and uncover the tree adaptation strategies during different phases of prior-wildfire tree growth and ecosystem recovery. Post-fire changes in the climatic response of tree rings have to be considered when reconstructing climate, as trees damaged by wildfire have to be excluded from the analysis, which is a common practice in dendroclimatology. On the other hand, well-defined wildfire-induced signatures in multiple tree-ring parameters could help in detecting and reconstructing fire history even for trees, stands and ecosystems with no clear detection of fire damage, like fire-scars.
Our reconstruction of the post-wildfire history for the forest ecosystem on permafrost testifies a multidecadal effect of forest fires on ASL thickness, ground vegetation and tree growth. This effect is well pronounced in variability of all the tree-ring parameters, including stable isotope ratios. The δ 13 C and δ 18 O data were also shown to provide an added value to understanding the physiological response of trees to the wildfire impact in southern Europe (Battipaglia et al 2014, Niccoli et al 2019, or even demonstrate environmental and physiological consequences of fire deficits for trees in the Western United States (Voelker et al 2019). Our study, for the first time in the boreal zone, demonstrates opportunities to unravelling wildfireclimate-permafrost-tree interactions when using treering isotopes. Similar studies need to be expanded to different tree species across the permafrost zone of boreal forests.