- Research Paper
- Published:

# A new approach for predicting board MOE from increment cores

*Annals of Forest Science*
**volume 78**, Article number: 78 (2021)

## Abstract

###
*Key message*

**Increment cores can provide improved predictive capabilities of the modulus of elasticity (MOE) of sawn boards. Multiple increment cores collected at different heights in a tree provide marginally increased accuracy over a single breast-height core, with higher labour costs. Approximately 50% of the variability of the static bending MOE of individual boards is explained by the predicted MOE obtained from a single increment core taken at breast height.**

###
*Context*

Prediction of individual board MOE can lead to accurate optimisation of the value extracted from forest resources, and enhanced decision-making on the management and allocation of the resource to different processors, and improve the processors ability to optimise grade allocation.

###
*Aims*

The objective of this study is to predict the MOE of individual sawn boards from the MOE measured from cores collected from standing trees.

###
*Methods*

A five-parameter logistic (5PL) function and radial basis function interpolants are used to obtain a continuous distribution of MOE throughout a log. By developing a “virtual sawing” methodology, we predict the individual board MOE for sixty-eight trees consisting of locally developed F_{1} and F_{2} hybrid pines (*Pinus* *caribaea* var. *hondurensis* × *Pinus**elliottii* var. *elliottii*).

###
*Results*

Moderate correlations for individual board predictions are observed, with *R*^{2} values ranging from 0.47 to 0.53. Good correlations between average predicted board MOE and average measured MOE are also observed, with *R*^{2} ≈ 0.83. A pseudo-three-dimensional approach, accounting for variation in height in the tree, affords marginally greater accuracy and predictive capability at the cost of increased data collection and processing. By using a single breast-height core, we can obtain a similar level of prediction of individual board MOE.

###
*Conclusion*

We have presented a novel non-destructive evaluation approach to predict the MOE of individual boards sawn from trees. This approach can be adapted to other wood properties, and other wood products obtained from trees.

## Introduction

Most commercial softwood shows a large variation of wood properties radially (from the pith towards the bark) which are generally more pronounced than variation within the same growth ring along the stem (Zobel and Van Buijtenen 1989). In softwoods, growth rings near the pith usually consist of a large proportion of earlywood tracheids (which have larger diameter and thinner cell walls than latewood tracheids), that gradually transitions to a larger proportion of latewood tracheids with increasing ring number from the pith (Liu et al. 2019).

The wood located near the pith, which generally shows relatively large gradients in properties, is called corewood or juvenile wood, whereas the wood located outside this zone is called outerwood or mature wood (Liu et al. 2019; Cown 1992; Zobel and Sprague 1998) . The corewood has shorter cells with smaller diameter and thinner cell walls, higher microfibril angle, lower specific gravity, lower modulus of elasticity (MOE) and strength compared to outerwood. Although corewood has some positive attributes for the paper and pulping industry, it presents significant problems for sawmills producing structural timber because of lower MOE, larger shrinkage and stability. In this work we focus on the commercial softwood species collectively known as the Queensland southern pines, that are commonly used for structural purposes.

The standard structural grade ranking of the Queensland southern pine resource is dictated by its stiffness as characterised by the measurement of the MOE (Australian and New Zealand Standard AS/NZS 4063.1:2010 2010b), since this resource is stiffness limited. Stiffness limited indicates that if the board stiffness achieves the required grade, the strength will exceed this grade, and hence, the MOE defines the mechanical grade of the board (or other structural products) extracted from the log. Non-destructive measurement tools that accurately measure board MOE allow enhanced genetic selection, site matching, harvest planning schedules, improved allocation of the resource to different processors and facilitate improved processor settings, product performance and grade recovery.

There are various methods and tools to estimate log MOE; however, there has been limited work to predict individual board MOE using data that can be obtained via non-destructive methods. Our main hypothesis is that it is not possible to predict individual board MOE from log-level information such as average log MOE. Therefore, we require a method to predict board MOE from non-destructive evaluation techniques.

Increment borers are used to extract cores from living, dead and felled trees for analysis of growth trends of tree ring patterns and to evaluate wood quality. The coring method is relatively inexpensive, rapid and simple. Increment coring is the most widely used sampling technique for wood density analysis (Gao et al. 2017), and has been used to develop models of wood density across a range of species (Kimberley et al. 2015; Fries and Ericsson 2006; Pokharel et al. 2014; Jordan et al. 2008).

In addition to density modelling, increment cores have also been used to analyse other properties and factors affecting wood quality, such as MOE, fibre coarseness, fibre wall thickness, annual ring growth and microfibril angle (Hong et al. 2015), to help decision-making for the improvement of juvenile wood (Lenz et al. 2011), and how wind regime affects the incidence of compression wood (Dinulica et al. 2016). Wood density and MOE variation in boreal softwoods (black spruce, balsam fir, jack pine) and hardwoods (paper birch, trembling aspen) were estimated using near-infrared spectroscopy on 30,159 increment cores from 10,573 inventory plots (Giroud et al. 2017).

Ivković et al. (2008) collected 12-mm bark-to-bark increment cores at breast height (1.3 m) for basic density analysis and bark-to-pith cores were assessed by SilviScan to obtain individual ring value MOE, and ring-area weighted averages. The objective of their study was to examine the variability and relationship between stiffness, strength, shrinkage and basic wood properties. They found that the variability in wood stiffness and strength from pith-to-bark was very high with greatest change near the pith.

Increment coring has been widely used for density analysis and developing correlation between properties. However, limited studies have been conducted to predict wood value with a focus on MOE. For example, in southern pine, Harding (2008) used 12mm diameter increment cores to assess extracted and unextracted basic density, spiral grain and microfibril angle from a range of southern pine genetics and silvicultural trials. He reported moderate correlations between core MFA measured on Silviscan and stress-wave velocity of the Fakopp and Wood Spec instruments. Kain (2003) also used a 12-mm-diameter increment core to assess spiral grain and variation in earlywood and lower latewood density and the ability to select for density in southern pine trials.

The IML Resistograph PD400 (IML System GmbH, Wiesloch, Germany) has been successfully used to measure the basic density of timber (Downes et al. 2018). Recently, Sharapov et al. (2019) have demonstrated the ability to obtain moderate correlations with static bending MOE of small test samples of clear timber using the Resistograph; however, it is currently not sufficiently reliable for measurement of MOE, hence our preference for increment coring in this work.

We use cores extracted from transverse discs taken from destructively sampled trees to construct models describing the pith-to-bark variation of MOE in the Queensland southern pines (*Pinus elliottii* and *P. caribaea* and their hybrids). We note that these cores are not specifically increment cores; however, the models developed are independent of the method used to collect the data. Furthermore, this approach can be classed as a non-destructive evaluation (NDE) technique, as outlined by Schimleck et al. (2019). We then use these models to obtain predictions of the stiffness of the boards sawn from these trees. This work is part of a large-scale collaborative project aimed at characterising the variability in commercial Queensland southern pine plantations, located in southeast Queensland and northern New South Wales.

Our goal is to evaluate the capacity to predict, using mathematical modelling, an individual board MOE with a low quantity of data (e.g. 5 to 12 data points per radius), and to this end we seek to model both the pith-to-bark and longitudinal variation of MOE within a tree. A description of the data has been given in Kumar et al. (2021), and is also discussed in Section 2.1. The radial variation of wood properties is an ontogenetic characteristic and occurs within each tree. It is a consistent pattern, and hence can be modelled. Tangential variation generally occurs in response to some external influence on the tree, such as environmental and climatic events (Zobel and Van Buijtenen 1989). Therefore, it is not possible, given our goals, to model the tangential variation. Furthermore, for MOE this tangential variability will be substantially less than the radial variability; hence, here, we only consider radial and longitudinal variation.

This will allow us to obtain a prediction of the MOE at any given radial position and height within the tree, and hence obtain a prediction of MOE for a product (structural board in this work) sawn from the tree. We note that instead of radius from the pith we could consider cambial age as governing the property variation, as done in our previous work (see Baillères et al. (2019) for details). When we are considering individual trees, there is little difference in the results obtained between radial position and cambial age. However, cambial age would be required when comparing between trees.

We describe a pseudo-three-dimensional approach, where we calculate the board MOE from four increment cores taken at different heights within the tree. This approach requires data collected at different heights, allowing us to determine whether longitudinal variation of MOE is significant for predicting individual board performance. We will then consider a pseudo-two-dimensional approach that utilises a single core obtained near breast height and neglects longitudinal variation in the tree, and compare this with the full pseudo-three-dimensional approach.

The methods developed here can predict the individual board MOE that can be obtained from a tree without actually cutting a tree. This enables important interventions at a grower level, such as pre-allocation of timber products to better match primary processors needs. It allows growers to identify plots that will never produce high value products so that these can be diverted to an appropriate facility and the site replanted to produce higher value products. Conversely, plots producing high value products early can be harvested at a younger age resulting in increased rotations from the same piece of land. Overall, the methods may be helpful to improve overall return, planning of investment, silviculture, harvest, appropriate resource allocation and marketing activities.

The remainder of this article is structured as follows. In the next section we discuss the data collection process, and the data obtained. We then outline our approach for modelling the pith-to-bark and longitudinal variation of MOE within the Queensland southern pines, before discussing our “virtual sawing” approach used to extract boards from a virtual log. We then compare our predicted board MOE values with those obtained by the standard four-point bending method (AS/NZS 4063.1:2010 2010b).

## Material and methods

### Data collection

The data (Psaltis et al. 2020) collected for this work come from the destructive sampling of three plots: 30 trees from within a typical commercial plantation (388 stems per hectare), and 38 trees across two spacing trials (between 200 and 2660 stems per hectare), including a Nelder wheel design (Nelder 1962), that covers a large range of growing conditions (Baillères et al. 2019).

Within each plot, the range of size classes were stratified to capture information about the full range of logs available within the plots. For the commercial plot, the diameters vary from 28.2 to 44.8 cm, and their heights from 23.2 to 27.8 m. The spacing trial plots have a diameter range from 21.2 to 51.5 cm and heights from 20 to 28.5 m. This ensures robustness of the predictive capability of our model.

This work focuses on the three major Queensland southern pine taxa; PCH (P. *caribaea* var. *hondurensis*), PEE (P. *elliottii* var. *elliottii*) and hybrid pine (PEE × PCH, both F_{1} and F_{2} hybrids) across age ranges: thinning age (15 to 20 years old) and harvest age (25 to 36 years old).

From the destructively sampled trees, four bark-to-bark cores were collected from transverse discs taken at 0.92 m, 2.34 m, 6.46 m and 7.88 m, as shown in Fig 1. The longitudinal variation of wood properties is much less than the radial variation, particularly above 5 m (Megraw 1985); therefore, we believe that collecting more data longitudinally will not substantially improve the predictive capability of the model. Additionally, each tree was merchandised into two peeler billets (see Psaltis et al. (2018) for details) and a 3.9-m sawlog. The sawlogs were weighed and measured and green stiffness properties determined using BING (Beam Identification by Nondestructive Grading) (Paradis et al. 2017), a resonance acoustic tool. The sawlog was then sawn and the boards were dried in accordance with industry recommendations for nominal structural framing dimensions (96 × 40 mm and 72 × 40 mm), using industry-provided sawing patterns. Full-length boards were tested using BING to obtain MOE data and measured for distortion (twist, spring and bow), then sub-sampled for static bending using the standard four-point bending test (AS/NZS 4063.1:2010 2010b).

The ring age and locations for each of the increment cores (all trees) were measured and recorded, and the cores were then segmented into 20-mm sections. Each segment was measured for MOE, density and shrinkage. The MOE of the segments was measured using an ultrasound time-of-flight approach. A computer code was developed in MATLAB^{®;} to take scanned images of the tree ring markings, detect their locations, and then divide them into the appropriate segments. Each of the measured wood properties were then attached to their corresponding segment. It is this segment data that forms the basis of our models for the board stiffness. For more details on the data collection and measurement approaches, the reader is referred to Kumar et al. (2021).

### Fitting radial sigmoid to increment core data

To model the pith-to-bark variation we assume that the MOE value for any given radial position can be represented in a functional form. It is well recognised that the MOE within an individual tree shows significant pith-to-bark variation, and exhibits a logistic growth profile (Baillères et al. 2005; Filipescu et al. 2014; McGavin et al. 2015). Therefore, the form of the function we assume to model a property, *p*, in terms of radial position, *r*, is the five-parameter logistic function (5PL), given by Gottschalk and Dunn (2005) and Dunn and Wild (2013)

where *α*_{0} is the maximum value of the property as *r* gets large, *α*_{1} is the minimum asymptote, *α*_{2} is the inflection point in the function, *α*_{3} controls the rate of variation of the property with radius, and *α*_{4} is the asymmetry parameter. We have chosen the 5PL as it allows for asymmetry in the data, and can often more accurately represent biological data than a four-parameter logistic function (Gottschalk and Dunn 2005; Dunn and Wild 2013).

Our aim is to then calculate the parameters, ** α** = (

*α*

_{0},

*α*

_{1},

*α*

_{2},

*α*

_{3},

*α*

_{4}), that give the best fit to the measured data. To do this, we consider a least squares approach. If we assume that the only error between our model,

*p*(

*r*), and the data,

**y**= (

*y*

_{1},…,

*y*

_{n}), is due to measurement error, we have the ordinary least squares problem.

However, if we assume that there is also error in the measurement of the independent variable, *r*_{i}, we have

where \(\hat {\delta }_{i}\) is the true error in the independent variable. This approach is known as orthogonal distance regression (ODR) (Boggs and Rogers 1990). ODR is preferred over standard regression techniques when there may be error in both the measured property and error in the independent variable (Boggs and Rogers 1990). In this work, the independent variable is the radial position, and this is associated with the measured MOE of a given increment core segment which may have errors or inaccuracies associated with its measurement. We show in Kumar et al. (2021) that ODR always gives a lower root mean square error when fitting the 5PL to our segment data.

Given both errors in measurement, *𝜖*_{i} and *δ*_{i}, the orthogonal distance between the curve and the point (*x*_{i},*y*_{i}) is given by

This leads to the unconstrained minimisation problem

In this work, we introduce additional constraints on the parameters, ** α**, by enforcing bounds to ensure physically meaningful results, such as restricting the maximum and minimum asymptotes to be within the physically possible MOE values of timber (1 to 25 GPa). To solve Eq. (4) subject to these additional constraints, we have employed MATLAB

^{®;}’s constrained minimisation routine fmincon (MATLAB Optimization Toolbox 2016a). Given an identified set of parameters, the function

*p*(

*x*) represents the variation of MOE with radius.

### Radial basis function surface fitting

We seek a model for the MOE that is continuous in both radial position and height, namely *s*(**x**), where **x** = (*r*,*z*) is the position in the tree in terms of radius and height. We use an approach based on radial basis functions (RBFs) (Carr et al. 2001) to obtain a smoothly varying surface for the desired property. RBFs have been used in a wide range of fields, from computer graphics to medical imaging (Carr et al. 1997). They provide a functional description of the data, and hence, the surface can be evaluated at any location. They do not require that data points be located on a regular grid, allowing for increased flexibility.

Our aim is to build a function that describes the property variation, given by

where *p*_{1}(**x**) = 1, *p*_{2}(**x**) = *r*, *p*_{3}(**x**) = *z*, *ϕ* is the radial basis function, and **c** = (*c*_{1},*c*_{2},…,*c*_{N}) and **d** = (*d*_{1},*d*_{2},*d*_{3}) are the unknown coefficients (to be determined). For details, the interested reader is referred to ??.

To apply this RBF approach to our problem, we assume our data, *y*_{i}, are samples from our fitted logistic functions. From the four cores we obtain information at four distinct heights within the tree, giving a dataset that captures the variation in both *r* and *z* for each individual tree in our study.

Following this process for each tree in the destructive sampling trials we have a description of the MOE at any location within each tree (bounded by the locations of the cores). This allows us to perform “virtual sawing” (discussed in the following section) to extract the timber boards from the pseudo-three-dimensional virtual log, calculate their MOE, and compare this value to the MOE obtained from static bending.

### Virtual sawing of standard boards

In this section we discuss our approach for extracting boards sawn from individual trees, to allow us to compare the predicted MOE with the corresponding measured MOE. The approach is identical whether we are considering a pseudo-three-dimensional or pseudo-two-dimensional approach, with the exception of the final calculation of board MOE. Here, we consider the pseudo-three-dimensional approach, where for two dimensions the mesh shown in Fig. 2 is reduced to a cross-section in the *x* − *y* plane.

To virtually saw boards from the log, we first require a digital representation of how the boards were sawn in the destructive sampling. The boards were sawn according to industry-supplied sawing patterns (output from a sawing optimiser); however, these were sawn in a small-scale mill in a government facility and may not accurately follow what would be obtained in a commercial mill. To replicate the pattern that eventuated when the boards were sawn, we have used log-end templates that were glued to the log prior to sawing. The location of the board centres could then be read from the templates, and these were manually transferred to clean copies of the template. An example template with board locations marked is shown in Fig. 3a.

We have implemented an image processing algorithm to identify the locations of the board centres from the scanned template images, and obtain a representation for each board. This process requires information on the board centres, what boards are aligned with each other, their dimensions and orientation. Each group of boards that is classified as being aligned is allocated an index, with index 1 being the reference alignment, and index 0 indicating no alignment. A line of best fit is calculated through the centre of the first alignment group, and the centres of the boards are then adjusted orthogonally to this line so that they lie on the line, as shown in Fig. 4.

We then compute an objective function that is composed of a shift along the line, a change in the dimension of the boards and the difference between the current board spacing and the saw kerf width (assumed constant). This difference is weighted more heavily, to ensure that the required saw kerf width is maintained following minimisation of the objective function. Transformation matrices are then computed for the rotation and translation of the boards. Once the first (reference) alignment group is completed, each subsequent group is adjusted, with a line of best fit taken as parallel to the reference line. The same procedure is undertaken to adjust the board centres and dimensions, and rotation and translation transformation matrices for each board are calculated. The final group of boards are those that are not aligned with other boards. Note however that their orientation will still be such that their faces are parallel/perpendicular to the reference line. We locate the nearest board to the current board, and adjust the dimensions and centre of the board to maintain the saw kerf width.

Figure 3b shows the digitised board layout for tree one, where we have two different board sizes sawn from this log (96 × 40 mm boards, denoted by A, and 72 × 40 mm boards, denoted by B). In this case, we have two rows of aligned boards (boards A1 to A5 and boards B1 to B5), with board A6 being unaligned. Following completion of this process, we then generate a canonical three-dimensional (or two-dimensional) mesh (in the *x*,*y*,*z* (or *x*,*y*) directions of Fig. 2) that discretises the dimensions of the board, with faces aligned with coordinate axes and centred at the origin. By applying the transformation matrix for each individual board to this canonical mesh, we are able to obtain a discrete representation for the volume (or cross-sectional area) of each board that has been sawn from a given log. Using the coordinates of each node in the mesh, we can calculate a radial (\(r = \sqrt {x^{2} + y^{2}}\)) and height position that can be used with our fitted surface, and hence calculate an MOE value at each node in the mesh. Note that in 2D, the MOE at each node can be computed directly from the fitted sigmoid functions.

We assume that surrounding each node there is a control volume of volume Δ*V*_{i} (or area in 2D), which takes the MOE value of the node at its centre, \(\mathcal {E}_{i}\). This is shown in Fig. 2. The volume averaged (Whitaker 1986) value of the MOE for a given board, \(\bar {\mathcal {E}}\), is given by

where *V* is the board domain, \(\mathcal {E}\) is the MOE throughout the board, and *V*_{T} is the total board volume. Equation (6) can be approximated using the discrete mesh as

where *N* is the number of discrete nodes for a given board mesh. Figure 5 shows an example of the digitised three-dimensional boards, where the shading represents the calculated MOE throughout the board. Note that these boards represent the test sample used for static bending taken from either the butt or top end of the 3.9 m sawlog. We can clearly see that the boards near the outer edge (bark) of the log exhibit higher values of MOE than those close to the pith.

To be able to compare the measured MOE values with those predicted from the virtual boards, we use the MOE obtained from the Australian Standard static bending test (AS/NZS 4063.1:2010 2010b) to give an MGP (machine graded pine) grading (AS/NZS 4063.1:2010 2010a).

The apparent modulus of elasticity in bending, *E*, is determined from the measurement of the change in vertical displacement, Δ*e*, at mid-span using the equation

where Δ*F* is the change in applied load. The static bending is performed on a test sample that is shorter than the full board, with the length of this test sample being determined by the cross-sectional dimensions of the board. For 96 × 40 mm boards, the test sample is 2 m long, whilst for the 72 × 40 mm boards it is 1.5 m. Each of these test samples is taken from the butt or top end, with the end being determined randomly.

## Results

Figure 6 shows the 5PL fitted to the measured segment data taken at four different heights within a given tree (core 1 to 4 in increasing order of height in the tree). Note that we have assumed radial symmetry here, allowing us to use data collected from both sides of the pith to fit a single curve. The fitted logistic curve gives a good representation of the observed trends in the data, whilst acting to reduce the variability. During early stages of a tree’s growth, the timber produced is highly flexible meaning low in stiffness. As the tree grows, it becomes increasingly rigid and this is reflected in the higher stiffness timber observed towards the outer edge of the tree (near the bark).

Figure 7 shows an example RBF surface fitted to the 5PL curves for a given tree. In this case, *λ* ≈ 0 (obtained by GCV), indicating our surface passes through the fitted 5PL curves. The measured segment data points have been overlayed in the figure to show the fit quality. This fitted surface allows us to smoothly connect the sigmoids obtained at different heights in the tree, to get an estimate of the property at any radius and height.

We now compare the results obtained using our two approaches. The pseudo-three-dimensional model has the advantage of capturing variation in height within the tree; however, the utility of this approach is limited as collecting this level of data is both time- and cost-prohibitive, and is difficult to perform as a non-destructive technique. The pseudo-two-dimensional approach utilises a single core extracted from near breast height, so may represent a more feasible non-destructive evaluation approach for the timber industry.

Figure 8 compares our predicted MOE using both the pseudo-two-dimensional and pseudo-three-dimensional approaches with those measured via the static bending standard test, for a range of representative trees. We see that for these trees our predicted board MOE captures the overall trends and ranges of board MOE measured from the static bending test. The pseudo-three-dimensional approach tends to predict higher MOE values for individual boards, and exhibits a lower root-mean-square-error (RMSE) than the pseudo-two-dimensional approach (mean of 1.911 versus 2.129) across all boards, when evaluated per tree.

Figure 9 shows a comparison between the measured MOE and the predicted values obtained using the pseudo-three-dimensional approach with radial basis functions fitted to the 5PL functions. We see that in each case (Commercial, Spacing/Nelder wheel, all trees) we obtain moderately good correlations for individual board MOE with *R*^{2} values between 0.48 and 0.53.

Figure 10 shows the correlation between the reconstructed board MOE using the pseudo-two-dimensional approach and static bending MOE across the commercial plantation and spacing trials, and all trees combined. Similar to the pseudo-three-dimensional approach, we obtain moderate correlations between the predicted and measured individual board MOE, with *R*^{2} values between 0.47 and 0.52.

Figure 11 compares the mean predicted board MOE with the mean measured board MOE across each of the trials, for both the pseudo-two-dimensional and pseudo-three-dimensional approaches. This shows that the pseudo-three-dimensional approach always exhibits a higher level of correlation with the mean measured board MOE than the pseudo-two-dimensional approach.

To evaluate model performance, we consider leave-one-out cross-validation (LOOCV), where a linear model is fitted to a training set consisting of all trees with one removed, and tested on the removed tree. This is done a total of 68 times, corresponding to each tree being in the test set. We have used the modelr (Wickham 2019) package in R 3.6.1 (R Core Team 2019). Using this approach we obtain a RMSE of 2.009 (with intercept) and 2.024 (without intercept) for the pseudo-two-dimensional method, and 1.915 (with intercept) and 1.919 (without intercept) for the pseudo-three-dimensional method.

This demonstrates that the pseudo-three-dimensional approach is generally more accurate in predicting individual board MOE when compared to the pseudo-two-dimensional method; however, as noted previously, the pseudo-three-dimensional method has a greater overhead in terms of data collection to obtain information at different heights within the tree.

Table 1 shows a comparison between the number of boards of grade MGP10+ (with MOE > 10 GPa, suitable for structural purposes), measured using the static bending method, and predicted from our pseudo-two- and pseudo-three-dimensional methods. In terms of number of MGP10+ boards, we can see that the best prediction is given by the pseudo-three-dimensional approach. However, the pseudo-two-dimensional method is able to predict 84% of the number of structural boards compared to approximately 90% for the pseudo-three-dimensional approach.

## Discussion

Our approach differs to many in the literature in that it is able to provide moderate predictive capabilities of *individual* board MOE, and very good predictive capabilities of the average board MOE obtained from a log. In Fig. 8, our predicted board MOE agrees well with the measured MOE. The range of values closely resembles the range measured by static bending, and we recover the trend of the stiffer timber (higher MOE) being found towards the outer edge of the tree. This is a direct result of our use of the 5PL function to model MOE with radius, and this behaviour is well recognised in the literature (Filipescu et al. 2014; McGavin et al. 2015). We expect there to be higher variability in the measured board MOE compared to the predicted MOE, and this is reflected in Fig. 8. This is because the predictions are made based on measurements of small, clear segments of wood that do not account for knots, and the 5PL will smooth any variation exhibited in the segment data.

When we regress the measured board MOE against the predicted MOE, shown in Fig. 9, the pseudo-three-dimensional model exhibits very little bias as the slope of the regression line is near unity in all cases. Furthermore, the coefficients of determination (*R*^{2}) show moderate correlations, comparable to alternative approaches in the literature (Launay et al. 2000; Wagner et al. 2003; Ishiguri et al. 2008). Furthermore, when using our pseudo-two-dimensional approach, we obtain similar *R*^{2} values, as seen in Fig. 10. Remarkably, by using a single core obtained from near breast height, the correlation again displays only marginal bias for all but the commercial plot, i.e. the predicted MOE is directly proportional to the static bending MOE as the slope of the regression line is very close to unity in each case, thanks to the analytical approach. We note that the location of the core used for the pseudo-two-dimensional approach is approximately 1.4 m below the lowest point of the sawn boards, and this could be the cause for increased bias compared to the pseudo-three-dimensional approach.

There are limited examples in the literature of comparing individual board MOE with measured values. Therefore, we have also considered the predicted average board MOE from each of our two approaches, compared to the average measured board MOE. These are obtained by simply taking the arithmetic mean of the measured and predicted board MOE for each tree.

The average board MOE has been compared to various non-destructive evaluation methods by a number of authors. Launay et al. (2000) developed a technique to measure the longitudinal MOE of a tree trunk and found moderate correlations between the tree MOE and average board MOE (*r* = 0.54). However, this corresponds to an *R*^{2} value of approximately 0.29; hence, the methods developed here represent a substantial improvement.

Wagner et al. (2003) compared stress-wave velocity measurements from standing trees to the average board dynamic MOE from Douglas-fir trees. Taking longitudinal stress-wave velocity measurements at 1.2 m (approximately breast-height), they found low coefficients of determination (*R*^{2} between 0.2 and 0.393). They were able to increase the correlation by taking additional measurements at multiple heights within the tree, and by combining average transverse and longitudinal stress-wave velocities they showed an *R*^{2} value of 0.591. They then considered the seven trees with the least variation in dynamic MOE, and were able to obtain a high *R*^{2} value of 0.839. We have been able to achieve similar or better results by taking a single breast-height core and analysing the MOE from individual 20 mm segments. When we consider our pseudo-three-dimensional approach utilising measurements at four heights within a tree, we obtain a substantial improvement over the approach of Wagner et al. (2003) (*R*^{2} between 0.69 and 0.85).

Ishiguri et al. (2008) investigated the use of stress-wave velocity and Pilodyn to predict mechanical properties of boards sawn from Japanese larch trees. They showed strong correlations between the average board static bending MOE and stress-wave velocity of six trees (*R* = 0.834, *R*^{2} = 0.696), and moderate negative correlation between Pilodyn and the average board static bending MOE (*R* = − 0.563, *R*^{2} = 0.317). They then combined stress-wave velocity with Pilodyn measurements in multiple linear regression, to get a correlation coefficient of *R* = 0.864 (*R*^{2} = 0.746) for the average static bending MOE. This is similar to the results we were able to obtain from a single breast-height core.

One main drawback of the two approaches presented in this work is the number of cores that can be collected per day, and specifically for the pseudo-three-dimensional approach, the limited ability to incorporate cores taken higher from the tree into industrial processes. As outlined in Kumar et al. (2021), approximately 30 to 40 breast-height cores can be collected and processed per day. This is substantially less than the approximately 160 to 200 ST300 readings and 300 to 400 Resistograph traces that can be collected each day (Baillères et al. 2019; Schimleck et al. 2019); however, these techniques are unable to predict individual board MOE which is key to establishing the structural grade and optimisation of the sawing patterns to obtain enhanced value extraction. Furthermore, they are less accurate, and thus, the gain in productivity may be offset by the level of prediction when compared to our approach. Our method is transferrable to a range of species, as it requires no calibration to be performed (Rakotovololonalimanana et al. 2015). It requires limited equipment, and can be adapted to other measurement techniques that can obtain the radial variation of MOE within a tree, or other characteristics displaying continuous radial variation such as density, fibre length, and microfibril angle.

Our methods are able to obtain moderate correlations with individual board MOE. Nevertheless, there are a number of sources of potential error in our approach that can be improved on. A plausible explanation for the low *R*^{2} values observed in the study could be due to the low spatial resolution of the data sampling and that we ignored ring width, growth history and tangential variation in our models. Better spatial sampling may allow us to capture more of the variability observed in the measured data. Additionally, there are errors due to inaccurate location mapping of the sawing pattern of the boards, leading to incorrect values of board MOE being calculated. This is due to having to recover the locations after the boards had been sawn. Furthermore, there can be approximately 11% relative error in the measurement of the static bending MOE (Baillères et al. 2019), and the effect of this should be investigated. We also note that these models have been developed using data obtained from measurements of MOE on defect-free samples. The static bending MOE can be significantly affected by knots and defects present in the board, and accounting for this is an ongoing area of investigation.

Our results show that the pseudo-two-dimensional method provides a simple and accurate approach to predict the range of board MOE values obtainable from each log, and offers improved predictive capabilities over existing NDE techniques. Additional accuracy can be gained by collecting additional cores at varying heights within the tree; however, this may not be a substantial advantage compared to the additional data required.

## Conclusions

In this work, we have developed mathematical models describing wood property variation in the Queensland southern pines based on cores extracted from transverse discs. These methods are directly transferable to increment cores obtained from standing trees. We have focused on prediction of the modulus of elasticity (MOE) as the main wood property governing value. This is the primary limiting property for the structural use of this resource and the most important mechanical property for structural end-uses which has a direct impact on structural timber grade outturn and the value of the resource. The processes we have developed can be extended to other properties such as density, permeability, Modulus of Rupture (MOR, although this would be more challenging due to the localised nature of this property), provided reliable data are available.

We have shown that by utilising a single breast-height core we can obtain similar prediction capabilities for individual boards as using four cores from different heights within the tree. The pseudo-three-dimensional approach gives greater accuracy and higher *R*^{2} values when considering average board MOE, however has higher data collection costs.

## References

AS/NZS 4063.1:2010 (2010) Timber Structures Design methods. Standard, Australian Standard/New Zealand Standard

AS/NZS 4063.1:2010 (2010) Characterisation of structural timber part 1: test methods. Standard, Australian Standard/New Zealand Standard

Baillères H, Vitrac O, Ramananantoandro T (2005) Assessment of continuous distribution of wood properties from a low number of samples: Application to the variability of modulus of elasticity between trees and within a tree. Holzforschung 59:524–530

Baillères H, Lee D, Kumar C, Psaltis S, Hopewell G, Brancheriau L (2019) Improving returns from southern pine plantations through innovative resource characterisation. Tech. report, Forest & Wood Products Australia, Project number: PNC361-1415

Boggs P T, Rogers J E (1990) Orthogonal distance regression. Tech. rep., Departments of Commerce National Institute of Standards and Technology, Gaithersburg

Buhmann M D (2004) Radial basis functions: Theory and implementations. Cambridge University Press, The pitt building trumpington street, Cambridge

Carr J C, Fright W R, Beatson R K (1997) Surface interpolation with radial basis functions for medical imaging. IEEE Trans Med Imaging 16(1):96–107

Carr JC, Beatson RK, Cherrie JB, Mitchell TJ, Fright WR, McCallum BC, Evans TR (2001) Reconstruction and representation of 3D objects with radial basis functions. In: Proceedings of the 28th annual conference on Computer graphics and interactive techniques, pp 67–76

Cown D (1992) Corewood (juvenile wood) in Pinus radiata - should we be concerned? N Z J For Sci 22:87–95

Dinulica F, Marcu V, Borz S A, Vasilescu M M, Petritan I C (2016) Wind contribution to yearly silver fir (Abies alba Mill.), compression wood development in the Romanian Carpathians. iForest-Biogeosci Forestry 9(6):927

Downes G M, Lausberg M, Potts B M, Pilbeam D L, Bird M, Bradshaw B (2018) Application of the IML Resistograph to the infield assessment of basic density in plantation eucalypts. Aust For 81 (3):177–185. https://doi.org/10.1080/00049158.2018.1500676

Dunn J, Wild D (2013) Chapter 3.6 - calibration curve fitting. In: Wild D (ed) The Immunoassay Handbook (Fourth Edition), 4th edn. Elsevier, Oxford, pp 323–336. https://doi.org/10.1016/B978-0-08-097037-0.00022-1, https://www.sciencedirect.com/science/article/pii/B9780080970370000221

Filipescu C N, Lowell E C, Koppenaal R, Mitchell AK (2014) Modeling regional and climatic variation of wood density and ring width in intensively managed Douglas-fir. Canadian Journal of Forestry Research

Franke R (1982) Scattered data interpolation: Tests of some methods. Math Comput 38 (157):181–200

Fries A, Ericsson T (2006) Estimating genetic parameters for wood density of Scots pine (Pinus sylvestris L.) Silvae Genet 55(1-6):84–92

Gao S, Wang X, Wiemann M C, Brashaw B K, Ross R J, Wang L (2017) A critical analysis of methods for rapid and nondestructive determination of wood density in standing trees. Ann For Sci 74 (2):27. https://doi.org/10.1007/s13595-017-0623-4

Giroud G, Bégin J, Defo M, Ung CH (2017) Regional variation in wood density and modulus of elasticity of Quebec’s main boreal tree species. Forest Ecol Manag 400:289–299. https://doi.org/10.1016/j.foreco.2017.06.019, http://www.sciencedirect.com/science/article/pii/S0378112717306242

Gottschalk PG, Dunn JR (2005) The five-parameter logistic: A characterization and comparison with the four-parameter logistic. Anal Biochem 343(1):54–65. https://doi.org/10.1016/j.ab.2005.04.035

Harding K (2008) Resource Characterization of Slash Pine Plantations Wood Quality. Journal article, Forest & Wood Products Australia

Holmes C C, Mallick B K (1998) Bayesian radial basis functions of variable dimension. Neural Comput 10(5):1217–1233. https://doi.org/10.1162/089976698300017421

Hong Z, Fries A, Wu H X (2015) Age trend of heritability, genetic correlation, and efficiency of early selection for wood quality traits in Scots pine. Can J For Res 45(7):817–825. https://doi.org/10.1139/cjfr-2014-0465

Ishiguri F, Matsui R, Iizuka K, Yokota S, Yoshizawa N (2008) Prediction of the mechanical properties of lumber by stress-wave velocity and Pilodyn penetration of 36-year-old Japanese larch trees. Holz Roh Werkst 66(4):275–280. https://doi.org/10.1007/s00107-008-0251-7

Ivković M, Gapare W J, Abarquez A, Ilic J, Powell M B, Wu H X (2008) Prediction of wood stiffness, strength, and shrinkage in juvenile wood of radiata pine. Wood Sci Technol 43(3):237. https://doi.org/10.1007/s00226-008-0232-3

Jordan L, Clark Iii A, Schimleck L R, Hall D B, Daniels R F (2008) Regional variation in wood specific gravity of planted loblolly pine in the united states. Can J For Res 38(4):698–710. https://doi.org/10.1139/X07-158

Kain D P (2003) Genetic parameters and improvement strategies for the pinus elliottii var. elliottii × pinus caribaea var. hondurensis hybrid in Queensland Australia. PhD thesis, The Australian National University

Kimberley M O, Cown D J, McKinley R B, Moore J R, Dowling L J (2015) Modelling variation in wood density within and among trees in stands of New Zealand-grown radiata pine. N Z J For Sci 45 (1):22. https://doi.org/10.1186/s40490-015-0053-8

Kumar C, Psaltis S, Baillères H, Turner I, Brancheriau L, Hopewell G, Carr E J, Lee D, Farrell T (2021) Accurate estimation of log MOE from non-destructive standing tree measurements. Ann For Sci 78(8). https://doi.org/10.1007/s13595-021-01031-w

Launay J, Rozenberg P, Paques L, Dewitte JM (2000) A new experimental device for rapid measurement of the trunk equivalent modulus of elasticity on standing trees. Ann For Sci 57(4):361–359. https://doi.org/10.1051/forest:2000126

Lenz P, MacKay J, Rainville A, Cloutier A, Beaulieu J (2011) The influence of cambial age on breeding for wood properties in Picea glauca. Tree Gen Genomes 7(3):641–653. https://doi.org/10.1007/s11295-011-0364-8

Liu F, Xu P, Zhang H, Guan C, Feng D, Wang X (2019) Use of Time-of-Flight Ultrasound to Measure Wave Speed in Poplar Seedlings. Forests 10(8). https://doi.org/10.3390/f10080682, https://www.mdpi.com/1999-4907/10/8/682

MATLAB Optimization Toolbox (2016a) The MathWorks. Natick, MA

McGavin RL, Baillères H, Ferhmann J, Ozarska B (2015) Stiffness and density analysis of rotary veneer recovered from six species of Australia plantation hardwoods. BioResources 10(4):6395–6416

Megraw R A (1985) Wood Quality Factors in Loblolly Pine. TAPPI Press, Technology Park

Nelder JA (1962) New kinds of systematic designs for spacing experiments. Biometrics 18(3):283–307, http://www.jstor.org/stable/2527473

Paradis S, Brancheriau L, Baillères H (2017) Bing: Beam Identification by Non destructive Grading

Pokharel B, Dech J P, Groot A, Pitt D (2014) Ecosite-based predictive modeling of black spruce (Picea mariana) wood quality attributes in boreal Ontario. Can J For Res 44(5):465–475. https://doi.org/10.1139/cjfr-2013-0252

Psaltis S, Turner I, Carr E J, Farrell T, Hopewell G, Baillères H (2018) Three-dimensional virtual reconstruction of timber billets from rotary peeling. Comput Electron Agric 152:269–280

Psaltis S, Kumar C, Turner IW, Carr EJ, Farrell T, Brancheriau L, Bailleres H, Lee D (2020) Prediction of board moe from increment cores - data. University of the Sunshine Coast Research Bank [dataset], https://research.usc.edu.au/permalink/61USC_INST/1vg4fiv/alma99484008902621

R Core Team (2019) R: A Language and Environment for Statistical Computing. R Foundation for Statistical Computing, Vienna. https://www.R-project.org/

Rakotovololonalimanana H, Chaix G, Brancheriau L, Ramamonjisoa L, Ramananantoandro T, Thevenon M F (2015) A novel method to correct for wood MOE ultrasonics and NIRS measurements on increment cores in Liquidambar styraciflua L. Ann For Sci 72(6):753–761. https://doi.org/10.1007/s13595-015-0469-6

Schimleck L, Dahlen J, Apiolaza L A, Downes G, Emms G, Evans R, Moore J, Pâques L, den Bulcke J V, Wang X (2019) Non-destructive evaluation techniques and what they tell us about wood property variation. Forests 10

Sharapov E, Brischke C, Militz H, Toropov A (2019) Impact of drill bit feed rate and rotational frequency on the evaluation of wood properties by drilling resistance measurements. Int Wood Products J 10(4):128–138. https://doi.org/10.1080/20426445.2019.1688455

Wagner F G, Gorman T M, Shih-Yin W (2003) Assessment of intensive stress-wave scanning of Douglas-fir trees for predicting lumber MOE. For Prod J 53(3):36–39

Wahba G (1990) Spline Models for Observational Data. CBMS-NSF Regional Conference series in applied mathematics. SIAM

Whitaker S (1986) Transport in porous media, D. Reidel Publishing Company, chap Flow in Porous Media I,: A theoretical derivation of Darcy’s law, pp 3–25

Wickham H (2019) modelr: Modelling Functions that Work with the Pipe. https://CRAN.R-project.org/package=modelr, r package version 0.1.5

Zobel B J, Sprague J R (1998) Juvenile wood in forest trees. Springer, chap Characteristics of Juvenile Wood, pp 21–55

Zobel BJ, Van Buijtenen JP (1989) Wood variation: its causes and control. Springer Science & Business Media

## Acknowledgements

The authors also wish to thank the following for their valuable contributions to this work: Bruce Hogg, Tony Burridge and John Oostenbrink (Department of Agriculture and Fisheries, Gympie) for plantation sampling, log selection and harvesting, and Dan Field, Eric Littee, Rhianna Robinson, and Rica Minett (Department of Agriculture and Fisheries, Salisbury Research Facility) for sawing and testing.

## Funding

This work was supported by Forest and Wood Products Australia (FWPA) through project PNC361-1415, and partners University of the Sunshine Coast, Queensland Department of Agriculture and Fisheries, HQPlantations, Forestry Corporation of NSW, Hyne and Son Pty Ltd, Queensland University of Technology and Hancock Victoria Plantations. SP wishes to thank the Australian Research Council Centre of Excellence for Mathematical and Statistical Frontiers (project number CE140100049) for postdoctoral research funding.

## Author information

### Authors and Affiliations

### Corresponding author

## Ethics declarations

### Ethics approval

The authors declare that they follow the rules of good scientific practice.

###
**Consent for publication**

All authors gave their informed consent to this publication and its content.

###
**Conflict of interest**

The authors declare no competing interests.

## Additional information

**Handling Editor:** Jean-Michel Leban

### Availability of data and material (data transparency)

The datasets generated during and/or analysed during the current study are available from the University of the Sunshine Coast Research Bank, DOI: https://doi.org/10.25907/00007

### Code availability

The custom code and/or software application generated during and/or analysed during the current study are available from the corresponding author on reasonable request.

### Publisher’s note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Contribution of the co-authors** Conceptualisation: Henri Bailleres, David J. Lee. Methodology: Steven Psaltis, Chandan Kumar, Ian Turner, Elliot J. Carr, Troy Farrell, Loic Brancheriau, Henri Bailleres, David J. Lee. Formal analysis and investigation: Steven Psaltis, Chandan Kumar, Ian Turner, Elliot J. Carr, Troy Farrell, Loic Brancheriau, Henri Bailleres, David J. Lee. Writing— original draft preparation: Steven Psaltis, Chandan Kumar. Writing — review and editing: Ian Turner, Elliot J. Carr, Troy Farrell, Loic Brancheriau, Henri Bailleres, David J. Lee; Funding acquisition: David J. Lee. Supervision: Ian Turner, Troy Farrell, Henri Bailleres, David J. Lee.

## Appendix:

### Appendix:

Our aim is to build a function that describes the property variation, given by

where *p*_{1}(**x**) = 1, *p*_{2}(**x**) = *r*, *p*_{3}(**x**) = *z*, *ϕ* is the radial basis function, and **c** = (*c*_{1},*c*_{2},…,*c*_{N}) and **d** = (*d*_{1},*d*_{2},*d*_{3}) are the unknown coefficients (to be determined). There are numerous forms of radial basis functions available, such as linear, Gaussian, and multiquadric (Buhmann 2004). Here, we use a thin-plate spline radial basis function, given by Wahba (1990)

where *ρ* = ∥**x** −**x**_{i}∥ is the distance between a data point, **x**_{i}, (a centre of the RBF) and a point on the surface. Nonlocal bases, i.e. those where \(\phi (\rho ) \rightarrow \infty \) as \(\rho \rightarrow \infty \), may perform better than local bases. Furthermore, the thin-plate spline is not dependent on a user-set shape parameter (Holmes and Mallick 1998), and is invariant under translation and rotation transformations (Franke 1982).

To obtain the unknown coefficients, **c** and **d**, we begin by assuming the data, *y*_{i}, can be modelled as

where \({\epsilon }_{i} \sim {} N(0, \sigma ^{2})\) is the error. We must solve the minimisation problem,

where *λ* is the smoothing parameter and *J*(*s*) is the penalty functional, given by

For *λ* = 0, *s*(**x**) becomes a surface that interpolates the data, and as \(\lambda \rightarrow \infty \) we obtain the linear least squares solution (Wahba 1990).

Wahba (1990) shows that the solution to the minimisation problem (Eq. (A.4)) can be found by solving the linear system,

where *K* is the *N* × *N* matrix with *ij* th entry *ϕ*(∥**x**_{i} −**x**_{j}∥), *P* is the *N* × 3 matrix with (*i*,*k*) entry given by *p*_{k}(**x**_{i}), *I* is the *N* × *N* identity matrix, \({\square }^{T}\) is the transpose operator, **c** = (*c*_{1},…,*c*_{N})^{T}, **d** = (*d*_{1},*d*_{2},*d*_{3})^{T}, **y** = (*y*_{1},…,*y*_{N})^{T}, and **0** is the 3 × 1 zero vector.

We see from Eq. (A.6) that to account for the penalty term we simply adjust the diagonal elements of *K*. To compute the value of *λ*, we utilise generalised cross validation (GCV) (Wahba 1990). This involves minimising the GCV function, *V* (*λ*), where

*T**r*() is the trace operator, and *A*(*λ*) is known as the influence matrix which can be calculated from (Wahba 1990)

Here, *Q*_{2} is computed from the QR decomposition of *P*, namely

where \(Q_{1} \in \mathbb {R}^{N \times 3}\), \(Q_{2} \in \mathbb {R}^{N \times (N-3)}\) and *R* is an upper triangular matrix.

## Rights and permissions

## About this article

### Cite this article

Psaltis, S., Kumar, C., Turner, I. *et al.* A new approach for predicting board MOE from increment cores.
*Annals of Forest Science* **78, **78 (2021). https://doi.org/10.1007/s13595-021-01093-w

Received:

Accepted:

Published:

DOI: https://doi.org/10.1007/s13595-021-01093-w

### Keywords

- Non-destructive evaluation (NDE)
- Wood quality
- Modulus of elasticity
- Queensland southern pines
- Sawn boards