Comparison of models for estimating bark thickness of Picea abies in southwest Germany: the role of tree, stand, and environmental factors

Abstract

Key message

Bark thickness was shown to vary between regions, stands, and trees. Bark thickness prediction equations of different model complexity can be suggested depending on the purpose of application. Site and stand conditions, which influence variation of growth rate to a large extent, seem to have a strong influence on bark thickness, with better site quality leading to smaller relative bark thickness.

Context

For many applications in forestry and forest science, local or regional species-specific bark thickness equations are used to estimate inside-bark diameters with outside-bark diameter measurements.

Aims

The objectives of this study were (1) to assess variation in bark thickness due to tree and stand factors in two Norway spruce (Picea abies (L.) Karst) datasets from different time periods, (2) to compare and evaluate alternative established model forms for estimating bark thickness, and (3) to assess spatial variation in bark thickness to estimate the effects of environmental factors on bark thickness.

Methods

Different bark thickness models were chosen from the literature and compared for their predictive quality for new measurements and a dataset from the 1970s. Mixed-effect modelling was applied to account for the hierarchical data structure, and generalized additive mixed models were used to analyse spatial effects and the influence of climatic factors, such as precipitation and temperature.

Results

A strong positional autocorrelation of bark thickness within trees and within plots could be shown. Bark thickness was smaller in the new data compared to the measurements from the 1970s. The variation between stands could not be explained by the tested environmental factors, but tree age had a strong positive effect on bark thickness.

Conclusion

In the study region, the variation of site productivity and individual growth rate seem to have a strong influence on bark thickness, whereas no significant effect of large scale climatic factors could be found.

1 Introduction

Although tree and log diameters are usually measured outside bark, for many applications in forestry and forest science, diameters inside bark are more important. Diameters inside bark are used to calculate wood volume of logs or trees for trade and forest inventories, and also, bucking optimization software in harvesters depends on them. To convert outside-bark diameters to inside-bark diameters, bark thickness equations have been developed for many species in the last decades (e.g. Li and Weiskittel 2011). These can also be used to calculate bark volume, which increasingly gained interest with the shift in the commercial relevance of bark from an unwanted residue to a valuable fuel and a source for high-value biomaterials (Doruska et al. 2009). The prediction of available bark biomass is important to assess the potential of such technologies generating potential additional income for the forestry sector. Another important application for bark thickness estimates is growth and increment predictions that are based on outside-bark measurements.

A common definition of bark, which we also used in this study, includes all tissues outside the vascular cambium and comprises secondary phloem up to the last-formed periderm and the rhytidome, which comprises all layers of dead tissue outside of the currently active periderm (Martin and Christ 1970). Bark development depends on two lateral meristems: the vascular cambium, which forms the phloem, and the phellogen (cork cambium) in each periderm, which forms phelloderm to the inside and phellem (cork) to the outside (Evert 2009).

The first studies on bark thickness and bark volume in central Europe were performed in the late nineteenth century in Switzerland to develop conversion factors from log volume with bark to log volume without bark. This research by Flury (1897) was the beginning of still ongoing research on the development of bark equations for many species in many parts of the world. Bark thickness equations have been developed from simple linear regression models that predict the percentage of bark with only diameter measurements (e.g. Flury 1897) to more elaborate model forms that describe the bark thickness as a function of several predictor variables that are usually easy to asses on the stem. The coefficients of such equations are usually valid for one species for a certain region. For many coniferous species, it has been shown that bark thickness can be well described by tree-inherent variables, such as diameter outside bark, total tree height, and height of the measurement position (e.g. Li and Weiskittel 2011). For Norway spruce (Picea abies (L.) Karst), it has been shown that bark thickness increases with diameter outside bark and with relative tree height and is influenced by tree age, tree height, and tree form (Laasasenaho et al. 2005).

Numeral external factors that correlate with bark percentage of trees have been reported, and many of them can be interpreted in a way that slower tree growth leads to larger relative bark thickness, as the following examples illustrate. For Pinus sylvestris (L.), the latitude of forest stands, which mainly describes differences in temperature and the length of the growing season, could help explaining bark thickness variation (Wilhelmsson et al. 2002). A correlation between a low site index or yield class of stands and a larger relative bark thickness was shown for Norway spruce (Dimitrov 1976; Hoffmann 1958; Schmidt-Vogt 1986) and silver fir (Abies alba Mill.) (Božić et al. 2007). Sonmez et al. (2007) reported a positive effect of tree age on bark thickness for Picea orientalis (L.) Link. Contrasting results were reported for the broadleaved Nothofagus pumilio (Poepp. and Endl.) Krasser by Cellini et al. (2012), where relative bark thickness was lower on low-productivity sites.

Besides growth conditions, also, genetics can determine bark thickness. Provenance has been shown to influence the phenotypic development of bark thickness of Pinus contorta var. latifolia Engelm. (Persson and Downie 1992) and Pseudotsuga menziesii (Mirb.) Franco (Kohnle et al. 2012).

In this study, we focus on Norway spruce in central Europe. Measurements in the state of Baden-Württemberg, southwest Germany, in the 1970s revealed a spatial heterogeneity of spruce bark thickness in that region (Altherr et al. 1978). Silvicultural practices have changed in the study region within the last decades with trees planted at wider initial spacing and earlier thinning being performed. Additionally, more forests have been converted to uneven-aged or mixed-species stands. Changes in silviculture, changed climatic conditions, and increased nitrogen emissions altered tree growth trends within the last century in the study area (e.g. Yue et al. 2011) and potentially also changed the relation between diameter increment and bark development. Advanced technology for data processing and modelling is available today and allows us to compare the correlation structure within bark thickness data from the 1970s and the 2010s and analyse the new data for potential factors influencing bark thickness development.

The objectives of this study were to:

• 1. Assess variation in bark thickness due to tree and stand factors in two Norway spruce datasets from different time periods

• 2. Compare and evaluate alternative established model forms for estimating bark thickness

• 3. Assess the spatial variation in bark thickness to estimate the effects of environmental factors on bark thickness.

2 Material and methods

2.1 Data origin

Bark thickness data were obtained from two sources: an extensive dataset with 7712 trees from Altherr et al. (1978) and a smaller set with 508 trees from measurements acquired between 2014 and 2016. Dimensional details are listed in Table 1. Trees of both datasets were sampled throughout the 35,752-km² state of Baden-Württemberg in temperate forests with a broad range of site and stand conditions. For the new data, 508 trees were selected for measurements in 35 plots, covering the major growth regions within the state (Fig. 1). Elevation of the plots varied between 277 and 1219 m above sea level. Mean annual air temperature varied between 4.6 and 9.3 °C and annual precipitation sum between 806 and 1945 mm year⁻¹. Mixed-species stands and pure spruce stands were chosen, and trees from thinning operations as well as from final harvest in selection-cutting systems were chosen. Sample size for plots and trees per plot was chosen after the variation in the Altherr et al. (1978) data had been assessed by Stängle et al. (2016).

Table 1 Characteristics of the measured trees from datasets 1 and 2: diameter at breast height (dbh), double bark thickness at breast height (DBT_1.3), the total tree height, and age. Dataset 1 is from Altherr et al. (1978)

Fig. 1

Distribution of the sample plots of dataset 2. Sampling of dataset 1 was performed in the same area (plots not shown)

2.2 Measurements

The sampling for both datasets was performed as described by Altherr et al. (1974, 1978): trees were felled, delimbed, and measured in the forest before any further log manipulation was performed. Measurement locations were at breast height (1.3 m above ground) and along the tree bole in 2-m increments up to a top diameter of approximately 10–15 cm. This distribution of locations within trees was shown to be adequate to describe within-tree variability within the merchantable bole length above 1.3 m (Stängle et al. 2016). Variability of bark thickness below 1.3 m could, however, significantly increase due to changes in stem form and bark profiles. Diameter and bark thickness were measured twice (approximately perpendicularly) at each location using both a calliper and a Swedish bark gauge, respectively. Double bark thickness was calculated as the sum of the two bark measurements. Additional measurements for the new data were tree height, which was measured on the felled tree, and tree age. Tree height for the Altherr et al. (1978) data was modelled with mixed B-spline regression describing tree taper using the R package TapeR (Kublin and Breidenbach 2013).

2.3 Variability between stands and trees

Due to the nested sampling design, bark thickness measurements display similarity of measurements within each tree and within plots (e.g. Li and Weiskittel 2011). If the variance between groups (such as trees and plots in our case) is large, additional information on group-level predictors could explain this variance. The between-group variance can be expressed by the intraclass correlation coefficient (Schielzeth and Nakagawa 2013). Variance components of two linear bark thickness models of different complexity (see below) were analysed, and the intraclass correlation coefficients were calculated as the ratio of between-group variances to total residual variance. Depending on model form and the included covariates, intraclass correlation can differ substantially for the same data (e.g. Stängle et al. 2016). Equations (1) and (2) were chosen from the literature as both are established model forms for predicting bark thickness. Equation (1) was introduced by Loetsch et al. (1973) and was used by Zacco (1974) for Sweden and was suggested as useful for spruce species by Li and Weiskittel (2011, eqn 7). Equation (2) was introduced by Hannrup (2004) and is implemented in the harvester protocol Standard for Forest Data and Communication StanForD (Skogforsk 2012). Double bark thickness is predicted in the two linear models as follows:

$$ DBT={\beta_0+\beta}_1{d}_{\mathrm{ob}}+\varepsilon $$

(1)

$$ DBT={\beta_0+\beta}_1{dbh}_{\mathrm{ob}}+{\beta_2 d}_{\mathrm{ob}}+\varepsilon $$

(2)

where DBT is the double bark thickness at any point of the stem (mm), d _ob is the diameter outside bark at that position (mm), dbh _ob is the diameter outside bark at breast height (1.3 m height) (mm), ε is the residual of the model, and β ₀ , β ₁ , and β ₂ are the regression coefficients.

2.4 Model comparison

Six different model forms were compared for their fit quality and their predictive capacity for both datasets. To account for the hierarchical data structure, we applied mixed-effect modelling with random deviation of the intercept on the tree and the plot levels, leading to Eqs. (3), (4), (5), (6), (7), and (8):

$$ DBT={\beta_0+\beta}_1{d}_{\mathrm{ob}}+{b}_i+{b}_{i j}+{\varepsilon}_{i j k} $$

(3)

$$ DBT={\beta_0+\beta}_1{dbh}_{\mathrm{ob}}+{\beta_2 d}_{\mathrm{ob}}+{b}_i+{b}_{i j}+{\varepsilon}_{i j k} $$

(4)

$$ \log (DBT)={\beta_0+\beta}_1{d}_{\mathrm{ob}}+{\beta}_2 \log {d}_{\mathrm{ob}}+{\beta}_3{\mathrm{ind}}_{d bh}+{b}_i+{b}_{i j}+{\varepsilon}_{i j k} $$

(5)

$$ DBT={d}_{\mathrm{ob}}\left({\beta}_1+{\beta}_2\frac{h}{H}+{\beta}_3{\left(\frac{h}{H}\right)}^2+{\beta}_4 H\right)+{b}_i+{b}_{i j}+{\varepsilon}_{i j k} $$

(6)

$$ DBT={\beta_0+\beta}_1{d}_{\mathrm{ob}}+{{\beta_2 d}_{\mathrm{ob}}}^2+{{\beta_3 d}_{\mathrm{ob}}}^3+{b}_i+{b}_{i j}+{\varepsilon}_{i j k} $$

(7)

$$ \log \left( DBT/{d}_{\mathrm{ob}}\right)={\beta_0+\beta}_1{\left(1-\frac{h}{H}\right)}^{\beta_2}+{\beta}_3{\left(\frac{h}{H}\right)}^{\beta_4 H}+{\beta}_5{ d bh}_{\mathrm{ob}}+{\beta}_6\frac{H}{{ d bh}_{\mathrm{ob}}}+{b}_i+{b}_{i j}+{\varepsilon}_{i j k} $$

(8)

where ind _dbh is an indicator for breast height (if true = 1; else = 0), h is the height above ground at the measurement point along the bole (m), H is the total tree height (m), b _i is the random effect for the ith plot, b _ij is the random effect for the jth tree in plot i, ε _ijk is the residual error for the kth measurement in the jth tree in plot i, β _0 – 6 are the regression coefficients, and other variables are defined as above.

Equations (3) and (4) correspond to Eqs. (1) and (2), respectively, with the additional random effect terms. Equation (5) was found to be suitable for spruce by Wilhelmsson et al. (2002, eqn S5b), and Eq. (6) was introduced by Cao and Pepper (1986, eqn 4) and suggested by Li and Weiskittel (2011, eqn 4). Model forms of Eqs. (7) and (8) were described by Gordon (1983, eqs 9 and 3). The dependent variable of Eqs. (5), (6), and (7) in the original publications was the diameter inside bark and was changed for this study, so that predictions were made for double bark thickness for Eqs. (3), (4), (5), (6), and (7).

To account for positional autocorrelation within each tree, which was not eliminated by introducing the tree-level random effect, a first-order continuous autoregressive correlation structure (CAR1) was applied (Pinheiro and Bates 2000), which has widely been used in forestry (Weiskittel et al. 2011). This accounts for the similarity of bark thickness values within the tree along the stem, with closer positions being more similar to each other. A power variance function with d _ob as covariate was introduced to account for a larger residual spread for increasing diameters. Previous analyses showed correlated errors and heteroscedastic errors in the data and led us to introduce the above-described within-group correlation structure and within-group heteroscedasticity structure. Bayesian information criterion (BIC) was computed to quantify an improvement of the model fit by introducing those structures.

Fit statistics include mean absolute error (MAE) and root mean square error (RMSE). Fits of the models (6) and (8) were corrected after logarithmic back-transformation (Baskerville 1972; Sprugel 1983). For the evaluation of the predictive capacity of the models, tenfold block cross-validation was performed. Therefore, the data was randomly split in ten subsamples with each of them containing the same number of plots. In each of ten runs, one subsample was used as validation data and the remaining nine subsamples as training data. Results from all runs were averaged. Mixed model fitting was performed using the lme function of the package nlme (Pinheiro et al. 2016) in the software R version 3.2.2 (R Core Team 2015).

To check if the two datasets led to significantly different parameter estimates, the best-fitting model (8) was fit to both datasets separately. Residuals from both fits were compared with equivalence testing to see whether the model form could describe both datasets similarly well. Subsequently, both models were applied to calculate double bark thickness predictions for the new dataset and results were compared with equivalence testing. For both equivalence tests, the region of indifference was set to 0.5 mm, which defines the absolute size that the mean of the differences can reach so that H₀ is still rejected, and therefore, equivalence is assumed (Robinson and Froese 2004). The R package equivalence (Robinson 2016) was used to calculate a robust two one-sided t test (TOST) for equivalence of unpaired data for the residuals and for paired data for the bark thickness predictions, respectively.

2.5 Spatial variation of bark thickness

In the last step of this study, we tried to interpret the variation in bark thickness between stands, which could be caused by environmental factors. As exact stand position and tree age were only known for the newly collected data, this step could only be performed for dataset 2. To account for climatic effects, we included annual precipitation (mm year⁻¹) and mean annual air temperature (°C). The elevation of the sample plots (metres above sea level) and geographic coordinates (Gauss-Krüger projection, in m) were included to evaluate spatial effects that could be explained by climatic as well as by biotic or abiotic factors such as site index and soil materials, respectively. Additionally, we included tree age to account for growth rate of the single trees, which can be influenced by site quality and tree individual factors, such as suppression in their early growth phase by larger trees. Generalized additive mixed models (GAMMs) were chosen for this analysis, and predictor variables were introduced as thin-plate regression splines with the R package mgcv (Wood 2011). Geographic coordinates were introduced as tensor product splines to model a two-dimensional surface as suggested in Wood (2006). To account for positional autocorrelation within each tree, a first-order continuous autoregressive correlation structure (CAR1) was applied (Pinheiro and Bates 2000). The model had the following form:

$$ DBT=\beta +{f}_s\left({d}_{\mathrm{ob}}\right)+{f}_s\left({d bh}_{\mathrm{ob}}\right)+{f}_s(H)+{f}_s\left(\frac{h}{H}\right)+{f}_s\left(\mathrm{age}\right)+{f}_{te}\left( X\mathrm{coord}, Y\mathrm{coord}\right)+{f}_s\left(\mathrm{elev}.\right)+{f}_s\left(\mathrm{prec}.\right)+{f}_s\left(\mathrm{temp}.\right)+{b}_i+{b}_{i j}+{\varepsilon}_{i j k} $$

(9)

where f _s (d _ob), f _s (dbh _ob), f _s (H), ${f}_s\left(\frac{h}{H}\right)$, f _s (age), f _s (elev.), f _s (prec.), and f _s (temp.) are smooth functions of d _ob, dbh _ob, H, $\frac{h}{H}$, tree age (years), elevation, annual precipitation, and mean annual air temperature; f _te (Xcoord, Ycoord) is a tensor product spline of geographic coordinates; b _i is the random effect of the ith plot, b _ij is the random effect for jth tree in plot i, and ε _ijk is the residual error for the kth measurement in the jth tree in plot i.

Fit statistics of submodels of Eq. (9) were compared after dropping single explanatory variables using the dredge function of the MuMIn package (Barton 2016) to test significance of each model parameter. BIC was computed to rank the fits achieved by this procedure.

3 Results

3.1 Variability between stands and trees

For both datasets, a strong positional autocorrelation of bark thickness within trees and within plots could be shown. Total residual variance was similar between the two datasets (old and new data) and was lower for Eq. (2), which had one more predictor variable than Eq. (1) (Table 2). Differences in the correlation structure of the two datasets could be shown, as the largest proportion of variability was found between the plots in dataset 1 and between the trees in dataset 2. The intraclass correlation at the plot level in the old data was approximately twice as high as in the new data, which means that the variability within plots, and therefore between trees, was much lower in the old data. The within-tree variability was higher in the old dataset, which was indicated by a lower intraclass correlation coefficient at the tree level.

Table 2 Correlation structure of both datasets expressed by total residual variance and variances between plots, between trees, and within trees for the two tested models

3.2 Model evaluation

Fit statistics revealed that Eq. (5), which included d _ob and ind _dbh as predictor variables, performed best for model fitting, but generally worst in the validation process (Table 3). For prediction purposes, the nonlinear model form (Eq. (8)), which included ${d}_{\mathrm{ob}},{ d bh}_{\mathrm{ob}},\frac{h}{H},\mathrm{and}\ H$, was found to be the best model for both datasets. Equation (3) that only used the diameter outside bark as prediction variable was found to perform worse than Eq. (4), which included also the diameter at breast height. Parameter estimates of all models can be found in Table 4.

Table 3 Summary of fit statistics of mixed-effect models (3)–(8) and the respective prediction bias using only the fixed terms for validation. The best values per dataset are set in italics

Table 4 Parameter estimates with their standard error in parenthesis of the models (3) to (8) for both datasets

Fitting model (8) on both datasets separately led to equivalent prediction bias (p < 0.05, region of equivalence 0.5 mm), leading to the assumption that the model could describe both datasets equally well. However, predictions for the newly assessed data were significantly higher (p < 0.05, region of equivalence 0.5 mm), when the model had been fit on dataset 1 compared to predictions made with a model fit from dataset 2 (Fig. 2). Inverting the TOST revealed that a region of equivalence of at least 2.9 mm would have to be chosen, to assume equivalence of the bark thickness predictions made with equations that had been fit on the two datasets separately.

Fig. 2

Difference in double bark thickness predictions gained with two fits of Eq. (8): one fit on the Altherr et al. (1978) data (dataset 1) and one fit on the newly collected data of this study (dataset 2)

3.3 Spatial effects

The GAMM fitting of Eq. (9) and its submodels obtained the best fit (assessed with BIC) when only the tree-inherent predictors ${d}_{\mathrm{ob}},{ d bh}_{\mathrm{ob}},\frac{h}{H}, H$, and age were included. Geographic coordinates and the climatic factors could not significantly contribute in explaining the variation. Estimated effects of the tested explanatory variables on double bark thickness can be interpreted with their smooth curves, where each curve illustrates the partial effect for each variable on double bark thickness, which is the sum of all effects. The effect of diameter outside bark on double bark thickness increased almost in a linear pattern (Fig. 3a). The effect of diameter at breast height increased over the full range of measured trees (Fig. 3b). The effect of tree height decreased linearly (Fig. 3c). The relative height showed almost no effect in the region of 10% relative tree height; however, a strong positive effect below and a slightly increasing effect above 10% could be observed (Fig. 3d). The effect of age clearly showed thicker bark at higher tree age (Fig. 3e). Diameter at breast height naturally increases with age. The increasing effects of diameter at breast height and age could therefore, at least partially, be caused by collinearity. Figure 4, however, shows that both variables have positive effects on the bark thickness prediction at the same time, meaning that bark thickness in trees of the same dimension was higher when the trees were older.

Fig. 3

Estimated smooth functions describing the partial effect of the five covariates on double bark thickness (DBT): a) the smooth of diameter outside bark, (b) the smooth of dbh outside bark, (c) the smooth of tree height, (d) the smooth of relative tree height, and (e) the smooth of age. Shaded areas are 95%-confidence intervals

Fig. 4

Contour plot showing estimated double bark thickness in millimetres over age and breast height diameter outside bark (dbh_ob). Covariates that were included in the model but are not shown were set to the median of occurring values in the data. Shaded area approximately covers the region that was covered by observations

4 Discussion

The two analysed datasets showed clear differences, with dataset 1 from Altherr et al. (1978) having higher bark thickness in general and more variability between different plots and within trees. We attribute those differences in the correlation structure to a lack of standardisation of the measurement procedure for dataset 1, as the Swedish gauge is very sensitive to differences in pressure (Gray 1956; Mesavage 1969). For those measurements, different personnel was used for different plots, whereas the new measurements for dataset 2 were all made by the same operator. Other influencing factors might be sampling season, which was almost exclusively out of growth season for the new data and could not be reconstructed for the old data. Also, we tried to sample a possibly large range of tree diameters in each plot in the new data, which might have reduced similarity of measurements within plots and, therefore, reduced the part of the variation that could be explained by the nesting in plots.

Results showed a complex nonlinear model (Eq. (8)) to be the best model for bark thickness prediction. This model cannot be suggested for all applications, though. Depending on the purpose of application, different explanatory variables are available and different degrees of computational complexity are reasonable. If bark thickness is to be modelled for inventory data, for example, a complex model with many explanatory variables, such as Eq. (8), can be suggested because taper curves outside bark are often available. If total and relative tree heights are not known, as in the situation when a harvester measures a tree, the best possible model without these covariates should be selected (Eq. (3)).

Differences in bark thickness development between trees of different seed source origin have been shown by Kohnle et al. (2012) and Persson and Downie (1992), however could not be reflected in this study. Most Norway spruce stands in southwest Germany were planted, but the origin of the seeds or the reproductive material was often unknown or not documented (Hosius et al. 2006; Konnert 1991).

The bark thickness variation that was found between the plots in this study could not be linked to large-scale spatial patterns, such as climatic zones, because geographic coordinates did not explain the variation to a significant amount. Also, the climatic factors temperature and precipitation could not explain differences between plots. Age, however, was found to have a significant positive effect on bark thickness. These findings are in accordance with findings of Laasasenaho et al. (2005) and Sonmez et al. (2007). This leads to the conclusion that faster tree growth results in smaller relative bark thickness for Norway spruce. Increased site quality, which usually can be linked to faster growth of individual trees, was shown to reduce bark thickness also in other studies (Božić et al. 2007; Dimitrov 1976; Hoffmann 1958; Schmidt-Vogt 1986). Also, Laasasenaho et al. (2005) found that parameters, which described the individual growth, could explain differences between climatic regions.

It seems that no single environmental variable, such as climate or geography, can explain differences in bark thickness, but only the combination of many factors, which influence the productivity of a stand and the growth of individual trees. For bark thickness studies in the future, the assessment of tree age or site index is, therefore, strongly recommended to consider the effects of factors such as nutrient and water availability, aspect, or stand density.

The effect, that faster-grown trees have a thinner bark, might explain the thicker bark in dataset 1 compared to the new data. Since growth of Norway spruce in southwestern Germany has accelerated in the twentieth century (Yue et al. 2011), it can be expected that the trees of the new dataset had grown faster than the trees of dataset 1, which were harvested about 40 years earlier. The same equations with the same combination of independent variables showed to best explain bark thickness variation in the two datasets. This leads to the conclusion that the general processes influencing bark thickness were not different in two studied time periods.

Intraspecific bark thickness differences between sites could mainly have two physiological reasons: either differences in the cambial activity on the xylem and the phloem sides or differences in the activity ratio of cambium and phellogen. Phloem increment in Norway spruce seems to vary between sites (Gričar et al. 2015b) but is more stable and less subjected to fluctuations in environmental conditions than xylem increment (Gričar and Čufar 2008; Gričar et al. 2007; Prislan et al. 2013). However, phloem cells collapse within few years (Evert 2009) and total bark thickness is largely determined by the thickness of the rhytidome. Factors influencing rhytidome development are difficult to estimate as there is not yet much research done on the complex processes of cell collapse and simultaneous inflation of other cell types (Gričar et al. 2015a). We may speculate that the higher investment in bark of slower-growing trees is evolutionary adaptive, as slower-grown trees need more time to reach reproductive age and, therefore, have a higher risk of mechanical damage and exposure to radiation and extreme temperature or other factors that a thicker bark might protect them from.

5 Conclusion

With our analysis, we confirmed that bark thickness varies for Norway spruce between single trees and between groups of trees that were sampled in different sites. Models of different complexity can be suggested depending on the purpose of application. An increase in predictive accuracy can be obtained by using models with more predictor variables; however, this increases the costs for data acquisition. In our study, we could not gain convincing evidence for geospatial variation of bark thickness development within the study area. In the study region, the variation of site productivity and individual growth rate seem to have a strong influence on bark thickness, with better site quality leading to smaller relative bark thickness.