Sustainable forest management in the Central Western Carpathians (Goat Backs, NE Slovakia): the role of climate change
SUPPLEMENTARY MATERIAL
Tomáš Hlásny1,2 • Ivan Barka1,2• Ladislav Kulla1• Tomáš Bucha1• Róbert Sedmák3 • Jiří Trombik2
1National Forest Centre – Forest Research Institute Zvolen, Department of Forest and Landscape Ecology, T. G. Masaryka 22, 960 92 Zvolen, Slovak Republic
2Czech University of Life Sciences Prague, Faculty of Forestry and Wood Sciences, Department of Forest Protection and Entomology, Kamýcká 129, 165 21 Prague 6, Czech Republic
3Technical University in Zvolen, Faculty of Forestry and Wood Sciences, T. G. Masaryka 24, 96092 Zvolen, Slovak Republic
Tomáš Hlásny
Ivan Barka
Ladislav Kulla
Tomáš Bucha
Róbert Sedmák
Jiří Trombik
SUPPLEMENTARY MATERIAL AStudy region and applied forest management
The Goat Backs Mts. model region is located in northeastern Slovakia (Central Europe) in the mountain range of the Central Western Carpathians (Fig. A1). The region covers an area of 8,226 hectares with 62.4% forest cover and an altitudinal range of 650–1,555 m a.s.l. Air temperature during growing season (April–September) ranges from 12–15°C, and growing-season precipitation ranges from 380–510 mm. Cambisols and Podsols prevail, while Rendzinas occur on the calcareous bedrock at the highest elevations.
Fig. A1 Location of the Goat Backs, Mts. model region in Europe (A, B), and regional forests, main soil types, and elevation zones (C). Such parameters were used to identifyforest stands representative of the study region to be used for forest development simulations
Since 2003, the incidental felling related to wind and snow damage followed by bark beetle outbreaks has been annually approaching 100% of the total volume of harvested timber. Even-aged management with the rotation period of 95-140 years (depending on species composition and site) is applied. Tending is applied at a thicket stage to reduce stand density to the 90% of the full stocking (the full stocking for a given site is specified using growth and yield tables). Between 2 to 4 thinning operations of variable intensity are applied, depending on stand density and site index. Regeneration system is uniform shelterwood in stands with fir and/or beech admixtures. In such stands, two-phase shelterwood cut at a maximum area of 3.0 ha with 4 - 6 harvest cycles is applied. In spruce monocultures, a small-scale clearcutting system with 3 harvest cycles is applied, with a maximum clear-cut area of 3.0 ha.
Because of high amounts of sanitary felling induced by the wind and bark beetles, the forest standsare stabilized by support to species other than Norway spruce by means of artificial regeneration (see Supplementary material B for the used species rates in the artificial regeneration). The proportions of each species in the artificial regeneration are site-specific. Spruce trees are harvested preferentially, while non-spruce species (e.g. larch or pine) are preserved to support their regeneration for a variable period of time, and are cut thereafter. Protection from ungulates during regeneration is required to ensure the survival of non-spruce species. Because dead or damaged trees support the growth of bark beetle populations (e.g. Hlásny and Turčáni 2013), such trees are rapidly removed, and almost no dead wood accumulates in the forest.
SUPPLEMENTARY MATERIAL BCriteria on the selection of representative stand types to be used for forest development simulations
Tab. B1 Criteria used to select the representative stand types (RSTs) for forest development simulations. Each combination of site and stand criteria is represented by a single RST, and each RST is representative of a certain part of the study region (Fig. A1). Species proportions used in the artificial regeneration in each RST are also indicated. Species codes: S – spruce, F – fir, L – larch, B – beech, M – maple
Species composition / Elevation (m a.s.l.) / Soil type / RST ID / Area(%)* / Regeneration (%)Spruce / 600-800 / Cambisols / 1 / 0.38 / S 50, F 10, P 10, L10, B 10, M 10
801-1100 / Cambisols / 2 / 10.79 / S 30, F 20, L 10, B 30, M 10
Podsols / 3 / 0.87 / S 40, F 20, L 20, B 20
1101-1500 / Cambisols / 4 / 4.42 / S 40, F 20, L 20, B 20
Podsols / 5 / 4.92 / S 60, F 20, L 10, B 10
Rendzinas / 6 / 4.19 / S 30, F 10, P 20, L 10, B 30
Spruce with larch / 600-800 / Cambisols / 7 / 1.86 / S 50, F 10, P 10, L10, B 10, M 10
801-1100 / Cambisols / 8 / 15.84 / S 50, F 10, P 10, L10, B 10, M 10
Podsols / 9 / 2.53 / S 40, F 20, L 20, B 20
1101-1500 / Cambisols / 10 / 7.32 / S 60, F 20, L 10, B 10
Podsols / 11 / 1.96 / S 60, F 20, L 10, B 10
Rendzinas / 12 / 3.86 / S 20, F 20, B 30, M 30
Spruce with larch and pine / 600-800 / Cambisols / 13 / 1.21 / S 30, F 20, L 10, B 30, M 10
801-1100 / Cambisols / 14 / 7.20 / S 50, F 10, P 10, L10, B 10, M 10
Podsols / 15 / 0.20 / S 50, F 10, P 10, L10, B 10, M 10
Spruce with larch, beech and fir / 600-800 / Cambisols / 16 / 0.38 / S 50, F 10, P 10, L10, B 10, M 10
801-1100 / Cambisols / 17 / 5.69 / S 40, F 20, L 20, B 20
Podsols / 18 / 1.21 / S 60, F 20, L 10, B 10
1101-1500 / Cambisols / 19 / 1.12 / S 60, F 20, L 10, B 10
Rendzinas / 20 / 0.19 / S 20, F 20, B 30, M 30
Mixture of spruce, larch, pine, fir, beech and maple / 600-800 / Cambisols / 21 / 1.84 / S 30, F 20, L 10, B 30, M 10
801-1100 / Cambisols / 22 / 15.52 / S 50, F 10, P 10, L10, B 10, M 10
1101-1500 / Cambisols / 23 / 3.18 / S 30, F 20, L 10, B 30, M 10
Podsols / 24 / 2.20 / S 60, F 20, L 10, B 10
Rendzinas / 25 / 1.11 / S 20, F 20, B 30, M 30
*per cent of the study region
SUPPLEMENTARY MATERIAL CDescription of climate change scenarios and of the baseline climate used to drive the forest development simulations
Tab. C1 Future climate statistics refer to the period 2071–2100. Although the data refer to elevation 1,000 m a.s.l., the projected change in air temperature and precipitation is almost uniform across the region. Model coding (c1–c5) follows the magnitude of Δ(°C) in the period 2071–2100. Abbreviations: MAT – Mean annual air temperature; PVS – Precipitation during vegetation season; Δ – Change between a scenario and the baseline climate
Scenario ID / RCM / Driving GCM / MAT (°C) / MAT Δ(°C) / PVS (mm) / PVS Δ (%)baseline / — / — / 7.0 / 0.0 / 451.0 / 0.0
c1 / DMI-HIRHAM5 / BCM / 9.0 / 1.9 / 495.5 / 9.9
c2 / DMI-HIRHAM5 / ARPEGE / 9.2 / 2.2 / 294.2 / -34.8
c3 / ICTP-RegCM / ECHAM5-r3 / 9.8 / 2.8 / 432.0 / -4.2
c4 / SMHI-RCA / ECHAM5-r3 / 10.3 / 3.2 / 381.0 / -15.5
c5 / HC-HadRM3 / HadCM3Q16 / 11.9 / 4.9 / 351.2 / -22.1
— / Average c1-c5 / — / 9.5 / 3.0 / 400.8 / -13.4
— / St. dev. c1-c5 / — / 1.5 / 1.0 / 66.7 / 15.3
SUPPLEMENTARY MATERIAL DImplementation of tree species-specific mortality rates in the forest dynamics model Sibyla
The records on forest damages in Slovakia (Source: National Forest Centre, internal data) in terms of volume of salvage harvests were used to calculate the mortality rates specific to a tree species and age class. The parameterisation data cover the period 1998-2009, and are based on randomly selected forest management units representing 59% of the forested area (1.14 mil. ha) and 60% of the growing stocks (273.8 mil. m3) in Slovakia. Although particular damage agents are specified in the source data, we used the aggregated mortality rates. The total harvested volume of dead or damaged trees relative to the total volume of a given species and age class represents the mortality rate used to drive the forest development simulations. Hence, the presented mortality rates include both natural tree mortality and mortality induced by disturbances.
Although the source data relate to a stand scale, we implemented the mortality rates in the Sibyla model so as each tree and age class in a stand receives a specific mortality rate. Although such a use violates the natural stand dynamics, it provides a realistic output when the data are spatially aggregated or simulation plots are representative of larger areas, as is the case of the current study.
As can be seen in Fig. D1, the mortality rates of spruce, larch, beech and oak are increasing with age, while pine and fir show the highest mortality rates at early ages. The latter fact is related to the frequent damages by snow and ice to pine, and by game to fir; the competition also increases the mortality rates for young developmental stages. The decrease in mortality rates at ages above ca 120 years is due to the fact that oldest stands in the parameterisation data are locatedin protected and protective forests with lower rate of damage as compared with commercial forests. However, trees older than 120 years occur in the study region only marginally.
Fig. D1 Mortality rates of main tree species used in forest development simulations in the Goat Backs Mts. model region. The values represent the probability a stand is destroyed in a certain age (pane a) or up to a certain age (pane b). In the pane a, the left y-axis is for spruce mortality rates, while the right y-axis is for the remaining species.
SUPPLEMENTARY MATERIAL EDefinition of forest development indicators addressed in the current study
Biomass productionper species (m3 ha-1 yr-1)includes roots, trunk, stem, branches, and leaves in both dead and living trees, harvested or not.
Total harvested biomass volume(m3 ha-1 yr-1) contains all biomass extracted during silviculture and harvesting operations
Harvested timber volume (m3 ha-1 yr-1) contains that part of total harvested volume that is classified in the Slovak assortment tables as timber wood assortments (Petráš et al. 1996)
Fuel and pulp wood volume containsthat part of the total harvested volumethat is classified in the Slovak assortment tables as fuel and pulp wood assortments
Tree species diversity was evaluated using the true diversity index (Jost 2006), which is the exponential form of Shannon`s diversity index :
where is the number of species, the basal area of species (m2; in this study calculated using the basal area of tree species with diameter ≥ 5cm) and
Hence the true diversity index is defined as:
Tree size diversity index was evaluated based on post-hoc index presented by Staudhammer and LeMay 2001:
where and are numbers of DBH and height classes present in a stand, is the basal area of DBH or height class (m2), and is the basal area of a stand (m2).For DBH, we use 5 cm classes, and for height 2 m classes (Cordonnier et al. 2013).
Dead wood volume (m3ha-1) includes standing dead trees with DBH 5 cm and lying dead wood originating from trees with DBH 5 cm at any decomposition stage.
Number of large standing dead and living trees (n ha-1)containsall large dead and living trees with a DBH above 70 cm for conifers and above 50 cm for broadleaves. Annual probability of dead tree downing was specified in Cordonnier et al. (2013), and the variable was modified accordingly.
SUPPLEMENTARY MATERIAL F Climate change effect on driving variables in the Sibyla model in the Goat Backs Mts. model region
Fig. F1 Climate change induced changes in the suitability scores of environmental variables driving the forest development in the used forest dynamics model Sibyla. The scores are calculated by transformation of original variables using the tree-species specific functions. The difference in suitability score between the period 2071-2100 and the baseline climate corresponding with the period 1961-1990 are presented. The future climate is represented by the average of 5 climate change scenarios described in the Supplementary materialC. Red zones represent the decrease of respective suitability score and deterioration of growing conditions, while green zones represent an improvement in growing conditions.
Driving variable codes: NOx – atmospheric NOx (ppb), CO2 – atmospheric CO2 (ppm), SoilN – soil nutrients (0-1), Vegs – number of days > 10° C (days), Tampl – annual temperature amplitude (°C), T – mean temperature during vegetation period IV-IX (°C), SoilM – soil moisture (0-1), P – precipitation totals during vegetation period IV-IX (mm), Arid – de Martone aridity index (mm °C-1)
SUPPLEMENTARY MATERIAL GInter-climate change scenario variability of evaluated forest development indicators
Fig. G1Inter-climate change scenario variability of simulated total biomass production (BP). The per cent change between five climate change scenarios and the baseline climate is presented for the current management regime (right pane) and the no management regime(left pane). Description of the used climate change scenarios is in the Supplementary material C. Periods: 0 – 2010, 1 – 2011–2030, 2 – 2031–2070, 3 – 2071–2100
Fig. G2 Inter-climate change scenario variability of the true diversity index (D) and the tree size diversity index (H). The difference between five climate change scenarios and the baseline climate is presented for the current management regime (right pane) and the no management regime (left pane). Periods: 0 – 2010, 1 – 2011–2030, 2 – 2031–2070, 3 – 2071–2100
Fig. G3Inter-climate change scenario variability of the simulated number of large standing trees per hectare (LSTN) and dead wood abundance in m3 ha-1 (DWV). The per cent change between five climate change scenarios and the baseline climate is presented for the current management regime (right pane) and the no management regime (left pane).Periods: 0 – 2010, 1 – 2011–2030, 2 – 2031–2070, 3 – 2071–2100
The inter-climate change scenario variability of the total BP (i.e. that of all species in the region) under the CMscenario ranged from 0 to -35% of the baseline-climate simulations; the greatest decrease was with c5, the warmest scenario (Tab. C1). The variability with NM ranged from 0 to -20% of the baseline-climate simulations.The greatest decrease was with the scenario c5 and with c2, the driest scenario (Tab. C1).
The inter-climate change scenario variability of differences between the climate change- and baseline climate-driven simulations of both D and H was much higher with CM than with NM (Fig. G2). The changes ranged from a minor decrease to a significant increase in D, and from a minor increase to significant decrease in H. LSTN showed larger variability with NM than with CM, and range from 0 to -80% of the baseline climate simulations. The inter-climate change scenario variability of DWV was higher with CM (ca. -20 – 30 %) than with NM (ca. ±20 %).
SUPPLEMENTARY MATERIAL HRelative importance of main forest development drivers
Because forest development in the region is predominantly driven by forest disturbances and management, and potentially by climate change, we determined the proportion of variability in simulated values accounted for by these drivers. We used a factorial analysis of variance with higher-order interactions based on the generalized linear model (GLM) to evaluate the relative effect of these factors (for more details, see Hlásny et al. 2014). To facilitate such evaluation, the simulations were run not only with CM and NM but also under a disturbance-free scenario (results of these simulations are not presented in this study). Hence, the statistical design had the following structure: the Management factor contained the levels “no management” and “current management”, and the Disturbances factor contained the levels “disturbances-yes” and “disturbances-no”. The Climate factor contained levels “no-climate change” and “climate change”; the latter level contained the simulations run under the five climate change scenarios (Tab. C1). The analysis was conducted in Statistica v. 12 (StatSoft Inc.).
The quality of GLM fit to the evaluated indicators as indicated by the R2 values was variable, which limits some of our interpretations (Tab. H1). The results show that D is mainly driven by forest management (39% of the total variability), and this agrees with the significant difference in D simulated under the two management systems (Fig. 3 in the article); the remaining portion of the total variability was not captured by the GLM. The relative effects on H and DW cannot be evaluated because of low R2 values (0.01 and 0.18, respectively). DWV, however, was significantly affected (p < 0.001) by management and disturbances. LSTN was well fit by the model (R2 = 0.75), and management, disturbances, and their interaction accounted for 39, 22, and 14% of the total variability, respectively.
The model fit for TBP was also very good (R2 = 0.72), and management, disturbances, and their interaction accounted for 30, 26, and 15% of the total variability, respectively. For TVH (R2 = 0.36), management, disturbances, and their interaction accounted for 22, 7, and 8% of the total variability.
The relative effect of climate change was insignificant for all variables, even though the results are based on the period 2070-2100, when the divergence among the five climate-change scenarios used was high (Supplementary material C).
Tab. H1 The percentages of the total variability in stand development indicators attributed to the main forest development drivers. The data refer to the period 2071-2100. The underscored values are statistically significant at p<0.001. Abbreviations: Clim – climate change, Man – management, Dist – disturbances, × – interactions, H – tree size diversity index, D – true diversity index, DWV – dead wood volume, LSTV – large standing tree volume, BP – total biomass production, TVH – total harvested volume
Forest development indicators / Forest development driversR2 / Clim / Man / Dist / Clim x Man / Clim x Dist / Man × Dist / Clim × Man x Dist
D / 0.4 / 0 / 0.39 / 0 / 0 / 0 / 0 / 0
H / 0.01 / 0.01 / 0 / 0.01 / 0 / 0 / 0 / 0
DWV / 0.18 / 0 / 0.13 / 0.05 / 0 / 0 / 0 / 0
LSTV / 0.75 / 0 / 0.39 / 0.22 / 0 / 0 / 0.14 / 0
TBP / 0.72 / 0.01 / 0.3 / 0.26 / 0 / 0 / 0.15 / 0
TVH / 0.36 / 0 / 0.22 / 0.07 / 0 / 0 / 0.08 / 0
REFERENCES
Cordonnier T, Berger F, Elkin C, Lämas T, Martinez M (2013) Models and linker functions (indicators) for ecosystem services. Arange Deliverable D2.2. Project Report. FP7-289437-ARANGE. 94 pp
Hlásny T, Turčáni M (2013) Persisting bark beetle outbreak indicates the unsustainability of secondary Norway spruce forests: Case study from Central Europe. Ann For Sci 70(5):481–491
Hlásny T, Barcza Z, Barka I, Merganičová K, Sedmák R, Kern A, Pajtík J, Balázs B, Fabrika M, Churkina G (2014) Future carbon cycle in mountain spruce forests of Central Europe: Modelling framework and ecological inferences. For Ecol Manag 328:55–68Jost L (2006) Entropy and diversity. Oikos 113:363–374
Petráš R, Halaj J, Mecko J (1996) Sortimentačné rastové tabuľky drevín /Assortment growth tables of tree species/. Bratislava. Slovak Academic Press. 252 p
Staudhammer CL, LeMay VM (2001) Introduction and evaluation of possible indices of stand structural diversity. Can J For Res 31(7):11 05-1115