Models to predict the radial variation of stiffness, strength, and density in planted noble fir, Norway spruce, western hemlock, and western red cedar in Great Britain

Abstract

Key message

This study compares the measured radial variation in wood stiffness, strength, and density of noble fir, Norway spruce, western hemlock, and western red cedar by developing mixed-effects models for each property using age as the explanatory variable. These models could be used to simulate the effect of rotation length and species choice on sawn wood properties.

Context

Timber production in Great Britain relies primarily on Sitka spruce. The use of multiple species is desirable to mitigate against biotic and abiotic risks posed to a single species. When considering alternative species, quantifying and modeling radial variation in wood properties is important to determine the potential for sawn timber production at a given rotation length.

Aims

To build empirical models for the radial variation in wood properties that can account for species.

Methods

Clear-wood samples were produced along radial transects in trees from four conifer species: Abies procera Rehder, Picea abies (L.) Karst, Tsuga heterophylla (Raf.) Sarg., Thuja plicata Donn. ex D.Don. Modulus of Elasticity, Modulus of Rupture, and density were measured on each species according to established standards. Mixed-effects models were built using ring numbers from the pith and species as explanatory variables.

Results

The same model forms could be used across the four species. Nonlinear models were developed for the Modulus of Elasticity and density. For the Modulus of Rupture, a linear model was most appropriate. The effect of species in the models was significant.

Conclusion

At similar rotation lengths, noble fir, Norway spruce, and western hemlock can produce timber with comparable properties to Sitka spruce. Overall, western red cedar would have worse properties for structural use.

1 Introduction

A primary objective of many planted forests is the production of sawn timber for construction. The production of sawn timber is important for economic purposes, and its use in construction ensures that carbon is stored for a long term, which is an important objective of contemporary forestry. For sawn timber to be used in structural applications by construction, there are quality requirements that the material must meet. Three physical properties typically characterize wood material for use in structural applications: Modulus of Elasticity (MOE), Modulus of Rupture (MOR), and density (e.g., Ridley-Ellis et al. 2016). Some non-structural uses of timber in construction, such as facades, also require some knowledge of stiffness, strength, and density for design. In Europe, these properties are used to define grading categories called “strength classes” (CEN 2016), each one defined by characteristic values, namely mean MOE (MOE_m) and the lower 5th percentiles for MOR (MOR_k) and density (ρ_k). Thus, grading is based on the description of a population, not individual pieces. In order for a population to be graded to a certain strength class, the three properties must achieve the required characteristic values of that strength class. If the characteristic values are achieved for the entire population, this is a 100% yield. Lower yields can be obtained by removing the lower quality material until the characteristic values are achieved for the remaining population, the passing part constituting the percent yield. Low yields are undesirable as there are negative consequences on the economy and for production efficiency. It is therefore important for tree growers to produce material that is suitable for their intended market. For example, in Great Britain the majority of structural timber is produced from Sitka spruce (Picea sitchensis (Bong.) Carr). The species is typically grown for 35–45 years (Moore 2011), and sawn-timber is graded to the C16 strength class (bending: MOE_m ≥ 8 kN/mm², MOR_m ≥ 16 N/mm², ρ_k ≥ 310 kg/m³) with essentially 100% yields. MOR_k and ρ_k usually achieve values suitable for higher strength classes (Moore et al. 2013), which makes MOE the limiting property and therefore that of most interest to determine the quality of British-grown spruce. There is evidence that this is also true for other conifer species grown in Britain (Gil-Moreno et al. 2016; Ridley-Ellis et al. 2016), but these are not yet as well characterized as Sitka spruce. The example of Sitka spruce given here is a wood supply chain where the resource is understood and the rotation lengths are suited to producing a particular product. Increasing the rotation length would allow higher strength classes to be produced (Moore et al. 2012) and decreasing it would lead to lower yields. For resources that are not well characterized, there is a requirement to estimate their suitability for an end-use as structural timber and to consider the possible influence of rotation length.

The radial variation in wood properties is particularly important when considering the effects of rotation lengths (e.g., Bao et al. 2001; Moore et al. 2012). Due to its large influence on forest products, the radial variation in wood properties has been widely studied in material from planted forests (reviewed in Zobel and Sprague 1998 and Lachenbruch et al. 2011) for many individual species. In general, for softwoods from planted forests, MOE and MOR increase from the pith to the bark (Bendtsen and Senft 1986; Moore et al. 2012; McLean et al. 2015). Often, empirical models describing the radial trend in MOE tend to incorporate a nonlinear increase from pith to bark that reaches an asymptote in the outerwood (e.g., Leban and Haines 1999; Auty et al. 2016). Models describing radial trends in MOR (e.g., Vincent and Duchesne 2014; Auty et al. 2016) usually have a similar functional form to MOE. Density can either increase, as is typically the case in the many pine species like radiata (Pinus radiata D.Don), loblolly (Pinus taeda L.), and Scots pine (Pinus sylvestris L.) (e.g., Xu et al. 2004; Mora et al. 2007; Auty et al. 2014), or decrease, at least initially, as is the case in most other conifers (e.g., Debell et al. 2004; Gardiner et al. 2011; Xiang et al. 2014). In species where there is an initial decrease in wood density from pith to bark, model forms are commonly “combined” or use multi-part functions that permit an initial decrease followed by an increase to an asymptote (e.g., Xiang et al. 2014; Kimberley et al. 2017; Auty et al. 2018). Irrespective of the wood property, the species considered, or the environment in question, a common feature of all models is that they include ring number from the pith, sometimes termed “cambial age,” as an explanatory variable. Some models include additional terms such as the cellulose microfibril angle for predicting MOE (Auty et al. 2016) or ring width for predicting wood density (Gardiner et al. 2011). A disadvantage of these additional terms is that they will not be readily available to forest planners, especially for species that are relatively uncharacterized in the intended growing environment. While radial variation is an important component of the overall variation in wood properties, it is not the only source of variation, differences between trees and between sites also make sizeable contributions to this total variation. Moore et al. (2013) showed that while radial variation accounted for most of the variation in MOE and MOR, tree-to-tree variation on a site was also important, particularly in the case of wood density. Generally, the variation attributed to site-to-site variation was around 25%. A limitation of the study by Moore et al. (2013) was that it only considered one species, Sitka spruce, and extending that to multiple species is relevant at the national level. The current study aims to build models of radial variation that only use age and species, two parameters that can be considered by forest planners as explanatory variables. As the models will be used to help plan future forests with a more diverse patchwork of tree species, it is of real practical importance that the sources of variation include species choice, a key decision that is made at the design stage of planted forests. Further, the aim is to use mixed models because these can account for and quantify the variation that is expected between sites and between trees within a site, which should ultimately allow for realistic predictions of the population level metrics, especially where the objective is not simply to predict the mean.

Quantifying and modeling the variation of wood properties, often for tree breeding or silvicultural studies, are usually performed using clearwood specimens. Clearwood specimens are more cost-effective to produce, easier to handle and test, and offer a higher temporal resolution with respect to annual rings than structural-sized timber. Further, clearwood specimens enable the study of the variability of wood properties with age without the complex influence of knots and other defects, which is often desirable for research purposes. However, the absence of knots presents difficulties in estimating the MOR of structural-sized timber from clearwood. This is because for structural-sized timber, failure (rupture) in bending will most frequently be due to the presence of knots (Rochester 1938 in Zhou and Smith 1991). Further, the size, nature, and frequency of the knots are correlated to the force required to break the timber piece (Divos and Tanaka 1997). Consequentially clearwood MOR can be, on average, around twice the MOR of structural timber (Zhou and Smith 1991; Lavers 2002). On the other hand, smaller differences are expected between the MOE in clearwood and timber, because knots have a relatively small influence on MOE (Samson and Blanchet 1992; Kretschmann 2010; Vikram et al. 2011). Finally, when measuring the density of structural-sized timber for grading a section free of defects (i.e. clearwood) is used to measure density; therefore, density is almost equivalent between the two specimen geometries. The material that is used in the current study has already been characterized at the population level as structural timber (Gil-Moreno et al. 2016a). Therefore, in the current study, this will permit an opportunity to evaluate the models that will be built against these population metrics. This is an important step in order to consider the real-world applicability of an approach that relies on clearwood specimens, and it will be informative for other researchers using the clearwood approach.

The objective of this study is to build empirical models of MOE, MOR, and wood density that are capable of quantifying the main variation in those properties in a way that is useful to forest planners aiming to introduce new species into the timber supply chain. The case study will consider four conifer species that currently form a minor component of the forest resource in the British Isles, but which are being considered for wider planting to increase the species diversity, and hence disease resilience, in the timber supply chain. The selection of species for the current study was primarily directed by their potential for their wider spread use in planted forests on the basis that some stands had already been grown to maturity, although they were not widely studied. The species are Noble fir (Abies procera Rehder), Norway spruce (Picea abies (L.) H. Karst), western hemlock (Tsuga heterophylla (Raf.) Sarg.), and western red cedar (Thuja plicata Donn ex D. Don.). Of these four species, Norway spruce is already part of the existing commercial spruce mixture, where it is mixed with Sitka spruce, although it has not been investigated independently in this environment. Limited data (Lavers 2002) suggest that noble fir and western hemlock should produce timber that has similar properties to spruce, while western red cedar, familiar on the import market due to its natural durability, would produce timber that has lower stiffness and strength. The current research on any of these species in this environment does not include information on how these properties change with age, which is required for informed decisions on the rotation lengths that are likely to be required to obtain material that can be used in structural applications. Our hypothesis is therefore that we will observe differences between the species considered, but that the general radial trend in wood properties will be comparable between species. We will test this by considering the significance of the inclusion of species in mixed models. Our second hypothesis is that we will be able to produce realistic population-level metrics of structural timber MOE and wood density using age and species as the explanatory variables, but that MOR will be overestimated due to the lack of knots in clearwood. We will test this by comparing our model estimates to the published values at a comparable age. The specific aims of the study are therefore: (i) to examine the variance of MOE, MOR, and density across four conifer species (ii) to model the variation of wood properties with cambial age, derived from ring number, using models that can be applied across all species (iii) to investigate the significance and effect of species using these models (iv) to test the models built for clearwood for prediction of structural properties on timber populations (v) to indicate the likely effect of rotation length on wood properties.

2 Material and methods

2.1 Materials

Noble fir, Norway spruce, western hemlock, and western red cedar were the species used in this study. The material and the study sites were previously described in detail in a study examining the grading potential for these species (Gil-Moreno et al. 2016a), and a simple table of tree characteristics (Table 1) and an overview are provided here. Each of the four species was sampled from even-aged single-species planted forests in three different growing regions in Great Britain: approximately in the south, middle, and north. This was intended to provide a sampling over a representative latitudinal range typical of where these species might be grown. At each site, three non-overlapping replicate plots were established covering the geographical spread of the stand and avoiding edge trees. Plots were at least 0.02 ha, and the plot diameter was incrementally increased until a minimum of 12 living trees of merchantable size were included within the boundary. Diameter at breast height (DBH) was measured on all living trees > 7 cm DBH, and the trees were ranked and grouped into quartiles by DBH. One random tree was felled from each of the three upper quartiles to represent the range of tree sizes likely to be processed for timber. This resulted in a total of nine trees per species felled at each site. An extra western hemlock tree was felled in the north site as a necessity of the felling operations, and therefore, it was included in the study. After felling, a log of approximately 5 m in length, with the bottom end at breast height, was cross-cut per tree. The upper 1.6 m of that log, which had its bottom end at approximately 4.7 m high in the stem, was designated for clearwood sampling. The lower part was used to produce structural-sized timber in the study by Gil-Moreno et al. (2016a). A 30-mm-thick central slab was cut from each of the 1.6-m-logs following a bark-to-bark diameter orientated to avoid defects. Ring numbers were marked on the top end of the slab to determine the cambial age of the clearwood specimens. The cambial age of a specimen is here defined as the mean ring number between the first and last ring from the pith contained within each specimen. From each slab, defect-free sections were identified in order to represent a complete radius from pith to bark. From the defect-free sections, specimens of approximately 22 × 22 × 300 mm (tangential × radial × longitudinal) were produced by sawing. Following air drying, the final sample target dimensions of 20 × 20 × 300 mm (tangential × radial × longitudinal) were obtained using a planer thicknesser. These sample dimensions correspond to BS 373 (BSI, 1957).

Table 1 Summary characteristics of the sampled trees

2.2 Mechanical testing

The clearwood specimens were conditioned to a constant mass in a controlled atmosphere of 65% RH and 20 °C. This corresponds to a nominal moisture content of 12%. The specimens were tested in a three-point bending configuration according to the BS 373 (BSI, 1957) “2 cm” specimen using a H5KT universal testing machine (Tinius Olsen LDT, Redhill, UK). Accordingly, the test span was 280 mm, and the speed was 6.6 mm/min. The orientation of the annual rings was parallel to the direction of loading. Prior to testing, the density was calculated from the mass and dimensions of the specimens. Following testing, the mass was measured before and after oven drying to finally determine the moisture content. MOE and density were adjusted to 12% MC following the adjustments given in EN384 (CEN 2018).

2.3 Statistical analyses

All statistical analyses were conducted within the R software for statistical computing (version 4.0.2, R Core Team 2020). Mixed model analyses made use of the nlme package (version 3.1–148, Pinheiro et al. 2020).

2.3.1 Variance components

The data included a hierarchical structure of tree within-plot within-site within-species. Variance components for each of the tree properties were calculated based on a simple random effects model (Pinheiro and Bates 2000):

$${y}_{ijklm}=\mu + {\beta }_{i }+ {\beta }_{ij }{+ \beta }_{ijk }+{\beta }_{ijkl }+{\epsilon }_{ijklm}$$

(1)

where y_ijklm is the observation of a wood property (MOE, MOR, or density) of the clear wood specimen m from the lth tree from the kth plot at the jth site of the ith species, µ is the overall mean, β_i, β_ij, β_ijk, and β_ijkl are the random effects of species, site, plot, and tree, respectively, and ϵ_ijklm is the residual variance attributable to the within-tree variation. Model residuals were checked for normality.

2.3.2 Radial models of wood properties

An aim of our analysis was to describe the radial variation in the measured wood properties using functions of cambial age across multiple species growing in a similar planted forest environment. These functions were parameterized to the data using a likelihood-based mixed-modeling approach (Pinheiro and Bates 2000) in which it was possible to both account for the hierarchical structure of the data (here tree, within-plot, within-site) and estimate the variance associated with each stratum.

In the models that follow, all α parameters are fixed effects of the functional forms. An interaction with species is considered with each fixed effect. Random effects are denoted as β parameters with character suffixes that designate the appropriate group as previously specified for the variance components analysis. Numerical suffixes are provided for convenience of identification in the tables of results. Random effects were allowed to act as additive terms for all of the models, meaning that they did not affect the parameters associated with the rate. In general, attempts to assign random effects to the rate parameters were met with non-convergence. Models were chosen from the literature, and an initial evaluation of some different functional forms was carried out in Gil-Moreno 2018. The models presented here offered equivalent or better-fit statistics (R ², AIC, RMSE) to those that were selected from multiple candidates.

For modeling the radial trend in MOE the model used was an exponential function of cambial age with intercept, previously used as a base model by Auty et al. (2016) for Scots pine:

$${MOE}_{ijklm}={\alpha }_{1, i}.{e}^{{\alpha }_{2, i}.{Age}_{ijklm}}+ {\alpha }_{3, i}+ {\beta }_{1,j}+ {\beta }_{2,jk}+ {\beta }_{3,jkl}+ {\epsilon }_{ijklm}$$

(2)

where age is the cambial age of the specimen, α₁ is the difference between the maximum value and the baseline, α₂ is the rate of change, and α₃ is the maximum value or asymptote; β_1,j, β_{2, jk}, and β_{3, jkl} are the random effects of site, plot, and tree, respectively, and ϵ_ijklm is the residual variance. An interaction with species i was considered with each fixed effect. It follows that the starting value (i.e., the stiffness of the wood in the center of the tree at ages approaching zero) can be obtained as the sum of α₃ and α₁.

The model used for MOR was a linear function of cambial age. Exponential models were investigated (Auty et al. 2016; Leban and Haines 1999) but a linear model gave the best fit:

$${MOR}_{ijklm}={\alpha }_{4, i}.{Age}_{ijklm}+ {\alpha }_{5,i}+ {\beta }_{4,j}+ {\beta }_{5,jk}+ {\beta }_{6,jkl}+ {\epsilon }_{ijklm}$$

(3)

where α₄ is the slope in a random slope model for the effect of age, and α₅ is the intercept; β_4,j, β_{5, jk}, and β_{6, jkl} are the random effects of site, plot, and tree, respectively, and ϵ_ijklm is the residual variance. An interaction with species i was considered with each fixed effect.

The model used for wood density was a simplification of the model proposed by Xiang et al. (2014) for black spruce (Picea mariana (Mill.) B.S.P.) and used by Auty et al. (2018) for Sitka spruce. This function couples an exponential term, to allow for the initial decrease in wood density, with the Michaelis–Menten equation to allow for the subsequent increase in wood density to an asymptote.

$${\rho }_{ijklm}={\alpha }_{6,i} \times {e}^{({\alpha }_{7,i}.{Age}_{ijklm})}+ \frac{\left({\alpha }_{8,i}.+ {\beta }_{7,j}+ {\beta }_{8,k}+ {\beta }_{9,l}\right). {Age}_{ijklm}}{{\alpha }_{9,i}+{Age}_{ijklm}}+ {\epsilon }_{ijklm}$$

(4)

where α6 is the initial value of wood density (ρ) and when the cambial age approaches zero, α₇ is the rate parameter of the early decrease in ρ, α₈ is the maximum or asymptote value in the outerwood, and α₉ is the rate of increase between the minimum value and this asymptote. β_7,j, β_{8, jk}, and β_{9, jkl} are the random effects of site, plot, and tree, respectively, and ϵ_ijklm is the residual variance. An interaction with species i was considered with each fixed effect.

The significance of species on the parameters of the three models was considered by comparing models with and without the inclusion of the fixed effect of species in a likelihood-ratio test (α = 0.05). Models were therefore fitted using maximum likelihood as this allows for the comparison of models with different fixed effects (Pinheiro and Bates 2000). The 95% confidence intervals of the estimated fixed effect parameters were used to assign the parameters for each species to groups and assess the differences.

Model performance was assessed based on both the explained variance, quantified using the coefficient of determination or R², and prediction error, quantified using the using the root mean squared error (RMSE), based on the linear relationship between observed and predicted values so that:

$${R}^{2}=1-\frac{\sum_{m=1}{({y}_{m}-{\widehat{y}}_{m})}^{2}}{\sum_{m=1}{({y}_{m}-{\overline{y} }_{m})}^{2}}$$

(5)

$$RMSE=\sqrt{\frac{{\sum }_{m=1}^{N}{\left({y}_{m}-{\widehat{y}}_{m}\right)}^{2}}{N}}$$

(6)

where y is the property and N is the number of observations. All model residuals were checked for normality.

2.3.3 Prediction of structural-sized timber from the clearwood models

In order to simulate populations of structural timber and compare these with published values obtained from the same trees in a previous study (Gil-Moreno et al. 2016a), we used the models to predict wood properties for each species. For this purpose, we required to create the age input for the radial models developed in the current study that was representative of the material in the study by Gil-Moreno et al. (2016a). This was done by drawing random samples from a probability density function of sample age. First, we created a normal (Gaussian) probability density function of age by species using the mean cambial age and associated standard deviation for each species (as reported in Table of Gil-Moreno et al. 2016a) Then, we randomly drew 1000 age values from these distributions for each species and used these as inputs to the radial models.

2.3.4 Simulation of the effect of rotation length on timber properties

For the final stage of our analysis, the models were used for a simulation of the mean MOE and wood density for different rotation lengths. We did not include MOR as we were not satisfied that the predictions would be realistic. To provide the age input for the models in this simulation, we considered the age input as 1,2…n, where n was the length of rotation in years (we used an arbitrary 20, 30, 40, and 50 years for illustrative purposes). The predicted mean value at each ring was obtained using the fixed effects only. The predicted mean value for each year was weighted by the relative cross-sectional for each corresponding ring and then summed to give the average at age n. For this purpose, rings were considered as concentric circles and their respective area was calculated accordingly. Radii was defined as the ring number. For each ring, the outer radius was considered as the ring number and the inner radius was considered as the previous ring number. The relative area is then the ratio of the ring area to the area of a circle with radius n. Actual ring widths or the ability to predict them were not available for these species and the use of this relative weighting was to ensure that the lower area of inner rings was not over-represented in the average.

3 Results

3.1 Summary statistics and variance components

The mean cambial ages of the specimens ranged from 1 to 39 in noble fir; 2 to 51 in Norway spruce; 2 to 60 in western red cedar, and 2 to 57 in western hemlock. Table 2 shows the mean values, standard deviations, and 5th percentile values for MOE, MOR, wood density, and cambial age for each species. Differences between the mean cambial age of the species are apparent and should be considered in any interpretation of the wood properties values. The differences observed in Table 2 are discussed in the following sections modeling the wood properties, which will allow the comparison of species over time. In general, Pearson’s correlation between MOE and MOR was very strong, and density was more strongly correlated to MOR than to MOE (Table 3).

Table 2 Summary statistics for cambial age, MOE, MOR, and wood density determined on clearwood specimens

Table 3 Pearson’s correlation (r) for the wood properties

The variance components for each property are given in Table 4. The largest component of the variance in MOE and MOR was within the tree. For wood density, the within-tree variation was large, but there was an equivalent component associated with variance between-species. The variation within-site (i.e., the sum of the variation attributed to differences between trees within a plot within a site) was generally greater than the variance between sites within a species, and summed approximately 75%, 65%, and 55% of the variation in MOE, MOR, and density, respectively. The variation attributed to differences between plots was negligible.

Table 4 Variance components (%) for MOE, MOR, and wood density clearwood specimens

3.2 Radial trends in wood properties and models

The three wood properties are plotted by cambial age for each species in Figs. 1, 2, and 3. For each property, the radial trend in each species is quite similar. In the case of western red cedar, the two mechanical properties do not increase further after the first 20–25 years, whereas for the other species, an asymptote has not been reached within the age range of the data studied. The radial trend in MOE can be described as nonlinear, the trend in MOR is approximately linear, and for density, the trend for each species is an initial decrease followed by an increase. The parameters of the three models for the radial variation in MOE, MOR, and density are given in Tables 5, 6 and 7, respectively.

Fig. 1

MOE plotted by ring number from the pith. Left, observed values. Different sites are represented with different plotting symbols (+ = north latitude, Δ = middle latitude, ○ = south latitude) and locally weighted regressions are plotted as trendlines. Right, marginal predictions from the model

Fig. 2

MOR plotted by ring number from the pith. Left, observed values. Different sites are represented with different plotting symbols (+ = north latitude, Δ = middle latitude, ○ = south latitude) and locally weighted regressions are plotted as trendlines. Right, marginal predictions from the model

Fig. 3

Wood density plotted by ring number from the pith. Left, observed values. Different sites are represented with different plotting symbols (+ = north latitude, Δ = middle latitude, ○ = south latitude) and locally weighted regressions are plotted as trendlines. Right, marginal predictions from the model

Table 5 Parameters obtained for the nonlinear mixed model for the radial variation in MOE

Table 6 Parameters obtained for the linear mixed model for the radial variation in MOR

Table 7 Parameters obtained for the nonlinear mixed model for the radial variation in wood density

3.2.1 MOE model

MOE is plotted by species and cambial age in Fig. 1 alongside lines representing the marginal predictions from the model. The parameters for the models are given in Table 5. The model with only fixed effects explained 69% of the variation in MOE, while an additional 16% was attributable to the random effects using the mixed model. Species was significant by inclusion in the mixed model (p < 0.001) when examined with a likelihood-ratio test. An examination of the model parameters showed that species generally differed from each other in at least one of the parameter estimates, with western red cedar most obviously differing from the other species and having lower MOE throughout the radial profile. Noble fir had a low initial value of MOE but increased to be almost equivalent to western hemlock in the outer part of the profile.

3.2.2 MOR model

MOR is plotted by species and cambial age in Fig. 2 alongside lines representing the marginal predictions from the model. The parameters for the models are given in Table 6. The model with only fixed effects explained 51% of the variation, while an additional 21% was attributable to the random effects using the mixed model. Despite some slight non-linearity near the pith, a linear model between MOR and cambial age was deemed appropriate over the range of the data. The inclusion of species and its interaction with age was significant (p < 0.001) in a likelihood-ratio test. Norway spruce and noble fir had a similar rate of increase in MOR with age, but Norway spruce had a significantly higher intercept term and therefore had consistently higher MOR across the range of cambial ages. Western hemlock had a significantly greater intercept term than the other three species, but a lower rate term than noble fir and Norway spruce. Western red cedar had an equivalent intercept term to Norway spruce, but almost half the rate of increase.

3.2.3 Density model

Density is plotted by species and cambial age in Fig. 3 alongside lines representing the marginal predictions from the model. The parameters for the models are given in Table 7. In terms of the model fits, using a mixed-effects modeling structure through the incorporation of random effects improved the explained variability. The model with only fixed effects explained 61% (RMSE = 34.4 kg m^-3) of the variation in density, while an additional 20% was attributable to the random effects using the mixed model. Species was significant by inclusion in the model (p < 0.001). According to the grouping that was determined on the 95% confidence intervals of the parameter estimates, all of the species had at least one parameter belonging to a different group, in other words, none was identical. Western hemlock had the highest density throughout the range of cambial age. It had a starting value (α6) and asymptote (α8) that were in the highest groups assigned based on the 95% confidence intervals. Noble fir and western red cedar had equivalent starting values to western hemlock, but they decreased more, with only noble fir returning to a value within the range of western hemlock. Norway spruce had the lowest starting value, but attained an asymptote that was not significantly different from western hemlock, that is in both species, it fell within the same group. Western red cedar had a significantly lower asymptote, assigned to a different group, than the other species and consequentially had lower density than the other species throughout the majority of the range of cambial ages.

3.2.4 Prediction of structural-sized timber population metrics from the clearwood models

Table 8 shows the values of the structural-sized timber published in Gil-Moreno et al. (2016a) compared to the predictions from the model built in the current study using clearwood from the same trees. Results show that, although the predictions are consistently lower than the observed values, the models made reasonable approximations of mean MOE, though in general, the predicted standard deviation was only about half of what was observed. Larger differences are found for mean MOR, where the predicted values were almost double what was observed. Additionally, there is a much lower standard deviation, which becomes a lot more pronounced if expressed as the coefficient of variation (not shown). This difference is translated to a large and unreasonable overprediction of the characteristic value, the fifth percentile, where the predicted values are around three times the observed values. For density prediction, the models offer very good predictions of the values in structural-sized timber with very similar means and comparable but slightly lower standard deviations, which corresponds to a small overprediction of the characteristic values.

Table 8 Observed population metrics for structural-sized timber published in Gil-Moreno et al. (2016a) compared to the predictions from the models developed in this study

3.3 Simulation of the effect of rotation length

Finally, the models were used to simulate the effect of rotation length on the mean MOE and density by species (Table 9). In general, longer rotations can be seen to offer an increase in the investigated wood properties, though this is less perceptible in the case of western red cedar. From the predictions, we can consider a comparison with the typical rotation of Sitka spruce, 40 years, and the MOE value of 8.00 kN mm⁻² required to obtain C16 grade timber. Norway spruce had the highest MOE of the four species, 8.69 kN mm⁻², and western hemlock was almost identical. Noble fir had achieved average MOE of 8.01 kN mm⁻² while western red cedar was substantially less at 7.04 kN mm⁻². The simulations predicted Norway spruce and western hemlock would have values in excess of 8.00 kN mm⁻² with a rotation length of 30 years. Using a rotation length of 50 years, Norway spruce and western hemlock would produce MOE_m values of C18 strength class. Mean wood density is not used in assigning strength classes (the 5th percentile is); nonetheless, the mean should be indicative of relative differences in the 5th percentile. For wood density, western hemlock ranked highest at any rotation length, followed by Norway spruce, then noble fir. Western red cedar had the lowest density of these species and also showed the least change with rotation length, though it did increase. Western hemlock had a lower relative change than Norway spruce and noble fir, but these three species showed a clearer increase in wood density with rotation length.

Table 9 Simulation of mean MOE and wood density for different rotation lengths in years

4 Discussion

In this study, the radial variation of MOE, MOR, and wood density was modeled with respect to cambial age (mean ring number from the pith) and species using data from physical and mechanical testing of clearwood specimens that followed a pith-to-bark profile. The models were used to test for significant differences between species and to simulate the population-level metrics that are used in the grading of structural timber. The predicted population metrics were compared to published population metrics of structural-sized timber obtained from the same trees. Finally, the models were used to simulate the effect of rotation lengths on the mean wood properties for each species in a way that is relevant to forest planning.

Our first aim was to examine the variance in each wood property by a variance components analysis. The hierarchy of the data included species at the highest level and within tree variation at the lowest level. For each property, results showed that even if examining different species, the variation within a single tree was the largest source of variation. We are not aware of other studies that have made this comparison across species. For Sitka spruce, Moore et al. (2013) also showed that within tree variation was the largest source of variation in a single species. In the current study, the proportion of variance attributable to within-tree variation was particularly high for MOE (71%). This is of special relevance because MOE_m is generally the characteristic value that limits the grading for timber in Great Britain. Thus, the results indicate that the length of time for which trees are grown has a bigger effect than the choice of species. For the other wood properties investigated, species had a higher influence, but it was only in density that the relative proportion was comparable to the within-tree variation. In the context of this study, the large proportion of variation that occurs within a tree justified the modeling approach used.

Our second aim was to model the radial variation in each of the wood properties. In the case of MOE, this followed a radial trend similar to that observed for other species such as Sitka spruce (McLean et al. 2015), Scots pine (Auty et al. 2016), loblolly pine and radiata pine, both reviewed in Schimleck et al. (2018). We are not aware of published radial models for any of the four species considered in the study. Some studies do exist with a focus on the radial variation of MOE in Norway spruce (e.g., Kliger et al. 1998; Koponen et al. 2005) or that describe some of the variation in western hemlock (Langum et al. 1995; Jozsa et al. 1998; Wang et al. 2001). In the current study, a single model form was successfully parameterized to the four species for MOE. This model has previously been used as a base model for Scots pine growing in the UK (Auty et al. 2016). It is not possible to directly compare parameters because the model in the earlier study was expanded with linear terms to include the effects of sampling height (a constant in this study) and an explanatory variable for microfibril angle (not measured here). Nonetheless, the model form is robust enough to be used for at least five species growing in Great Britain and is likely to be suitable for other conifer species too. In the case of MOR, we are not aware of published models for any of the species studied. Auty et al. (2016) used the same model for MOR as for MOE in their study on Scots pine, but a linear model appeared to be more appropriate to our data, which spanned a shorter time period. We believe it is not reasonable to expect the same linear trend indefinitely, and the addition of older material beyond the range of ages examined would likely describe a non-linear model reaching an asymptote. Thus, the extrapolation of these linear models to older rotation lengths beyond what is considered here is not recommended. The final property considered with respect to radial variations was wood density. This is also the most widely studied (reviewed in Saranpää 2003). There are examples in the literature of studies including radial variation for three of the species considered here (Norway spruce: Petty and Macmillan 1990; Jyske et al. 2008—western hemlock: Singleton et al. 2002; Debell et al. 2004—western red cedar: Cown and Bigwood 1971), including some models (Mäkinen et al. 2007; Moltenberg and Hoibo 2007). We found no published information on the radial variation of wood density in noble fir. The model applied was the same used by Auty et al. (2018) for Sitka spruce and Xiang et al. (2014) for black spruce. The parameterization of this model in other species in another situation demonstrates the robustness of this form. Overall, the fact that we observed comparable trends across species for each property and that we could use the same model form for each property proved part of our first hypothesis that the general radial trend would be comparable across species. In this instance, we believe this is because conifers are physiologically similar and it is expected that the physical and mechanical properties of trees from each species will all respond to a comparable growing environment.

The third aim was to investigate the differences in species using these models and the main part of our first hypothesis was that we would see differences between species. This turned out to be true and we did see differences between species. However, these were generally quite subtle. From the model fits, western red cedar had a lower MOE throughout the radial profile. Western hemlock and Norway spruce had similar radial profiles of MOE. Noble fir had lower overall MOE than western hemlock and Norway spruce, largely influenced by the lower values near the pith, almost half that in Norway spruce. These differences were reduced in the outer portion of the tree, and noble fir reached comparable stiffness to that in western hemlock. Western hemlock had higher MOR than the other species in the corewood, but Norway spruce and noble fir increased to similar values in the outerwood as it had a higher rate term in the model. MOR was lower in western red cedar than the other three species throughout the profile, indicating its wood may be less well suited to structural applications. Regarding density, western hemlock had the highest throughout the profile, followed by Norway spruce and noble fir, while western red cedar had the lowest values.

Our fourth aim was to use the models to predict population-level metrics for timber and to compare these with published values for the same material. In the case of MOE, the models slightly underpredicted the mean and underpredicted the standard deviation. One explanation for this could be that our observed population-level metrics for clearwood were lower than the associated population-level metrics observed for structural timber. We believe that the most likely explanation for this is the fact that different test methods were used, specifically structural timber is tested a four-point bending test with different bearing conditions. According to the study by Brancheriau et al. (2002), a larger shear deformation and indentation of the load point into the test piece in a three-point bending configuration will account for differences between results for the same specimens in three and four-point bending, with the results for three-point bending being lower by the same order of magnitude that we observed here. In the case of MOR, we did not expect that the prediction of structural timber MOR would work well, and it did not. It has been shown that the relationship between MOR in structural-sized timber and clearwood specimens directly produced from that timber is not good (Butler et al. 2016; Kranjnc et al. 2019). While this does not set a good foundation, our aim was not a direct prediction of specimen properties, but rather the population level metrics. In this respect, our models overpredicted mean MOR and underpredicted the coefficient of variation. While it might be possible to apply a heuristic correction to the mean and similarly heuristically expand the variance, our predictions would still not follow the same simple ranking as the observations if ordered by species. For the case study considered here, MOR is generally of less consequence when assigning a grade than MOE for timber produced from planted forests in Great Britain (Moore et al. 2013; Gil-Moreno et al. 2016a), there are examples where it can sometimes determine the timber grade (Gil-Moreno et al. 2019), and therefore, a solution is required. We propose that this solution will not come from the study of clearwood alone, but it is necessary to consider the relationship of MOE to MOR in structural-sized timber and if possible incorporate information about knots in structural timber (e.g., Divos and Tanaka 1997). It is possible to model knot dimensions with forest growth (e.g., Manso et al. 2020), and therefore, a future solution could be considered incorporating this. Finally, in the case of wood density, we expected this would work well and it did. Our models slightly underpredicted the coefficient of variation and therefore overestimated the indicated 5th percentile value. Estimating the extremities of a population is always likely to be difficult because it takes relatively very few values to influence the observed 5th percentile and relatively many samples to influence the mean. Thus, we feel it is more important that model predictions of the mean and coefficient of variation are reasonable, and in this case, the simple ranking of species was correct and we believe the models performed well. In summary, we can accept our hypothesis that we would be able to predict population-level metrics of MOE and wood density, but not MOR.

Our final aim was to simulate the influence of rotation length and species choice on wood properties. Given that our models did not perform satisfactorily for MOR, we presented this only for MOE and wood density. We chose to present mean wood density as it is more reliable than predicting the 5th percentile, but the 5th percentile could be estimated using the coefficients of variation if an overestimate could be accepted. Again, for the case study in question, this is the most relevant because in Great Britain the grade will be determined primarily by MOE. Based on the estimated mean MOE, western hemlock and Norway spruce should produce high yields of C16 with a rotation length of 30 years, while for noble fir a rotation length of 40 years is required. Western red cedar would produce lower yields of C16 than the other three species, but it can produce C16 timber. Further, these results quantify the investigated properties for use in non-structural applications such as facades. For these applications natural durability and visual appeal might be the primary reason for species selection; however, knowledge of the properties investigated in this study is nevertheless required.

These models and estimates can be combined with the appropriate growth and yield functions in order to help make decisions related to forest planning. In this context, it is important to consider the minimum rotation length that would be required to produce a particular grade as other factors such as age of maximum mean annual increment or wind risk will have at least an equivalent bearing on the decision about what to plant and where. Ultimately, it would be desirable to integrate the models produced here directly with growth and yield models to ensure that the area weighting of cross-sectional averages was correctly proportioned. We used a relative proportioning to demonstrate the principle and it is not possible to investigate if this is entirely adequate in the current absence of diameter growth functions for the species in question.

5 Conclusion

This study investigated the characterization of a forest resource for the three-key wood properties that characterise structural timber for grading using a study of clearwood specimens. Models for the properties of MOE, MOR, and density, using cambial age (tree ring number) as the sole predictor, that can be used across species were built for noble fir, Norway spruce, western hemlock, and western red cedar growing in planted forests in Great Britain. Despite some differences between species, all of the species examined can be used to produce structural timber, and possibly, the similarity in the underlying physiological processes is more notable than the differences in species per se. The cross-species applicability of the models suggests that these are very robust models for many conifer species growing in planted forests, and in the future, they could and should be expanded to include other species. The use of age as the explanatory variable allows these models to be combined with growth and yield simulators, should they become available for the species-environment combinations in question. This could enhance forest production forecasting in a way that is meaningful for wood processing. The current models could also be useful in economic studies that can attribute different values to different quality classes of wood for example. However, a weakness in the results that must be addressed with future research is that the model for MOR is unsatisfactory, both in the ability to predict the properties of structural-sized timber, and in that it was a linear model where a nonlinear trend is arguably more likely. To address this latter aspect, a further study would need to consider trees from substantially older forests, which do not currently exist in Great Britain. To address the issue about the prediction of the MOR of structural lumber, one path that could be investigated is to predict the knot characteristics that may occur with tree growth. Despite the shortcoming in the ability to reliably predict MOR from the models produced in this study, MOR is not a limiting factor in this use case and we were successfully able to use the models to demonstrate the influence of rotation on length on the wood properties, where longer rotations are shown to offer better wood properties.