Comparison of whole-tree wood property maps based on near-infrared spectroscopic calibrations utilizing data at different spatial resolutions

Sample origin

Two sets of 18 P. taeda trees, aged 13 and 22 years, respectively, were sampled from a half-sib progeny trial planted at the International Paper’s (IP) Southlands facility located near Bainbridge, GA, USA. The trees were a subsample of two larger populations selected to encompass the range of cellulose, lignin and specific gravity (SG) variation as measured by IP. Selected trees were felled and samples removed for wood property and NIR analyses. Four bolts (0.75 m long) per tree with each bolt cut in the center of each 25% of total height and discs (25 mm thick) were taken at 1.5-m intervals along the stem of each tree, giving 9–13 discs per tree (number depended on tree height). A total of 191 discs were collected from trees aged 13 years, while 220 were obtained from the 22-year-old trees, giving a total of 411 discs.

Sample preparation – radial strips

Pith-to-bark radial sections were cut from each disc using a bandsaw. Section dimensions were 12.5 mm longitudinally by 12.5 mm tangentially with the radial length corresponding to the radius of the disc. Radial sections were gently dried, glued into core holders and cut using a twin-blade saw to give strips 2 mm thick (tangentially) using the methodology described in Jordan et al. (2008).

Near-infrared spectroscopy

Schimleck et al. (2009) provided a detailed description of the methodology used to collect NIR spectra. Briefly, NIR spectra were collected in 10-mm radial increments from the radial-longitudinal face of each radial strip using a FOSS NIRSystems Inc. (Laurel, MD, USA) Model 5000 scanning spectrometer. The spectra were collected at 2-nm intervals over the wavelength range of 1100–2500 nm. The total number of radial spectra collected was 2569 (1114 from 13-year-old trees and 1455 from trees aged 22 years).

Wood property measurement

A subsample of 72 strips (two from each tree and representing different heights) was selected for SilviScan determination of ADD, MFA, MOE and several tracheid properties including coarseness (Evans 1994, 1999, 2006). From the same trees, four 0.75-m bolts were chipped using a 1.2-m Carthage chipper (Carthage Machine Company, Carthage, NY, USA) at the Institute of Paper Science and Technology (IPST), Atlanta, GA, and then thoroughly mixed to give a single composite sample per tree. Though several properties were measured for this study, our focus was on density and tracheid coarseness as they were the only properties measured on both radial strips and whole-tree samples.

The basic density (oven-dry weight, green volume) was measured on the whole-tree chip samples by water displacement and by oven-drying the chips (White et al. 2009). White et al. (2009) provided a detailed description of pulping. Briefly, for each whole-tree sample, duplicate cooks were performed on 1 kg of oven-dried whole-tree chips in a 10-l M/K digester. Two pulp grades (linerboard at 100 κ and bleachable at 30 κ) were obtained for each sample. For both pulp grades, the fiber (tracheid) coarseness was measured using a fiber quality analyzer (FQA) (OpTest Equipment Inc., Hawkesbury, Canada) (White et al. 2011); however, in this study, only the fiber coarseness for the linerboard pulps were used.

The measured wood properties (radial strips vs. whole-tree chips or pulp) were adjusted to a common basis to facilitate comparisons between them. The basic density values of the chips were converted to a SilviScan equivalent ADD which corresponds to a moisture content of 7.7%, and thus was a change from 0% to 7.7% moisture for the weight, and a change in the volumetric dimensions of 9%. The 9% change in dimensions corresponds to a change in volume from the fiber saturation point for loblolly pine, assumed to be 28.7% moisture content, to the 7.7% moisture corresponding to SilviScan ADD, using 12.3% as the total volumetric change for loblolly pine from the fiber saturation point (assumed 28.7%) to 0% for loblolly pine (Glass and Zelinka 2010).

Adjusting densities to a comparative basis is routine, but the same cannot be said for coarseness as during pulping, the tracheid width and wall thickness decrease due to the loss of hemicellulose and lignin (Scallan and Green 1975; Kibblewhite and Evans 2001). Reductions in tracheid dimensions change pulp coarseness and are related to reductions in pulp yield (Çöpür et al. 2005). To the authors’ knowledge, models for adjusting coarseness are not available in the literature for loblolly pine, so we used the radiata pine data presented in Kibblewhite and Evans (2001). Specifically, this paper does not present model parameters, but data for the individual trees are presented and thus utilized here to construct a linear model to adjust from coarseness measured on pulped fibers to coarseness measured on dry samples using SilviScan. The adjustment from pulped coarseness to dry coarseness is as follows:

where CoarsenessDry is the coarseness predicted on dry solid samples using SilviScan, and CoarsenessPulp is the coarseness measured on the pulp obtained from the kraft process (Kibblewhite and Evans 2001).

ADD and coarseness calibrations

Two distinct model types (spatial-specific and spatial-interpolated) were examined in this study. For the spatial-specific models (referred to in the following text as NIR-SilviScan models), SilviScan ADD and coarseness data were averaged over 10-mm sections from pith to bark for use in model building with the 10-mm NIR spectra data (477 in total, representing 179 and 298 spectra from trees aged 13- and 22-years, respectively) (Schimleck et al. 2018a). The 10-mm spatial-specific model was then used to predict ADD and coarseness for all 2569 spectra in 10-mm increments.

Spatial-interpolated models were developed at the whole-tree level. This required that the 10-mm spectral data from each tree be averaged to give a single representative spectrum (2569 10-mm increment spectra averaged to 36 in total) to allow for model building with the wood properties measured at the whole-tree level. Two approaches were used to determine the average tree spectra, with the first being an arithmetic average (NIR models based on these spectra are referred to in the text as NIR-Mean) for all spectra within a tree. The second was a basal-area weighted average (NIR models=NIR-Weighted) which weighed spectra based on the basal area they represented within the tree. For each wavelength from 1100 to 2500 nm in 2-nm intervals, the following steps were conducted. First, the basal area for each radial position represented within a disc was calculated as follows:

where area is the area that a specific 10-mm radial section represents within a disc, r_o is the outside radius of the 10-mm section, and r_i is the inside radius of the 10-mm section. Note that NIR spectra at the bark represent a greater area than NIR spectra at the pith, because their area is larger and thus they are given a larger weight within the radial strip. The weighted average spectrum for each disc was then calculated as follows:

where spectrum_disc is the weighted average spectra for a specific disc, area_i is the area that a specific 10-mm radial section represents within a disc, spectrum_i is a specific spectrum, area_disc is the total area of a disc, and n is the number of 10-mm radial sections within each disc. The weighted average spectrum for each tree was then calculated as follows:

where spectrum_tree is the weighted average spectrum for a specific tree, area_disci is the area that a specific disc represents, spectra_disci is the calculated spectra for each disc, area_tree is the total area for each tree, and k is the number of discs within each tree. Note that the weighted average NIR spectra at lower height levels within each tree (such as the stump) represent a greater amount of the total area of the discs than the NIR spectra at the tip of the tree.

Untreated spectral data, as used by Schimleck et al. (2009) and Mora and Schimleck (2009), were used to create the calibrations using partial least squares (PLS) regression. NIR-SilviScan calibrations were developed with four cross-validation segments while calibrations for the whole-tree data (NIR-Mean and NIR-Weighted average) were developed using leave-one-out cross-validation. The calibration performance was assessed using the following parameters:

1. Standard error of calibration (SEC) determined from the residuals of the final calibration;

2. Standard error of cross-validation (SECV) determined from the residuals of each cross-validation phase;

3. Coefficient of determination (R²), the proportion of variation in the calibration set that was explained by the calibration;

4. Ratio of performance to deviation (RPD_c) (Williams and Sobering 1993), calculated as the ratio of the standard deviation of the reference data to the SECV.

Determination of RPD allows comparison of calibrations for different properties that have differing data ranges and units; the higher the RPD, the more accurate the data is described by the calibration or predicted. The PLS models were then used to predict ADD and coarseness for all 2569 spectra in 10-mm increments. Figure 1 provides a summary of the three models.

Figure 1:

Summary of the three different models developed for air-dry density (ADD) and tracheid coarseness (C).

Maps of wood property variation within trees

Mora and Schimleck (2009) provided a detailed description of the methodologies they used to develop maps of within-tree variation for ADD and MFA. Three different methods were explored: Akima’s interpolation, Universal Kriging and semiparametric smoothing. Mora and Schimleck (2009) reported that for loblolly pine, maps generated by Akima’s interpolation method provided a good representation of the expected trends in ADD and MFA. In addition, they reported that a principal advantage of Akima’s algorithm over the Universal Kriging and semiparametric smoothing techniques was that only straightforward procedures were required, and there were no problems concerning computational stability or convergence. Owing to these features, only Akima’s interpolation method was used to provide the maps that were used to compare spatial-specific (NIR-SilviScan) and spatial-interpolated (NIR-Mean, NIR-Weighted) models. Maps for the 13- and 22-year-old trees were built using data predicted by the NIR-SilviScan, -Mean and -Weighted average approaches.

Comparison between predictions and measured data

For each age, probability density plots (where the area under the density curves has a probability of 1) were produced to show the variation in the predicted properties for the three different sets of ADD and coarseness values. The predicted NIR-SilviScan, NIR-Mean and NIR-Weighted average prediction values were compared for ADD and coarseness using a linear model to determine the model accuracy between the three methods. The whole-tree mean density and coarseness were calculated for each of the three methods by weighting the volume of the point compared to the overall volume. The weighted property was then summed up for each respective height and radial position combination to yield the mean ADD and coarseness of each tree.

Data analysis

All data were analyzed using the R statistical programing environment (R Core Team 2018) with the RStudio interface (RStudio 2018) and the tidyverse series of packages (Wickham and RStudio 2017). The NIR calibration models were developed using the pls (Mevik et al. 2013) and signal (Signal Developers 2013) packages. Mean and weighted spectra were calculated using the dplyr package (Wickham and Francois 2016). The ADD and coarseness within-tree variation maps (for the trees aged 13- and 22-years) were produced using the akima (Akima and Gebhardt 2016), fields (Nychka et al. 2015) and lattice (Sarkar 2008) packages.

Wood property calibrations

The NIR-SilviScan ADD calibration, developed using 10-mm SilviScan ADD data and NIR spectra measured in 10-mm increments, and the two spatial-interpolated ADD calibrations (NIR-Mean, NIR-Weighted), developed using whole-tree ADD and NIR data, are summarized in Table 1. The NIR-SilviScan ADD calibration was based on five factors and had a high R² (0.88) and good RPD_c (2.8). For the NIR-Mean and NIR-Weighted ADD calibrations, two factors were used and the models had very similar statistics. Calibration statistics were inferior to those reported for the NIR-SilviScan ADD model, but this could be expected as the range of the calibration data for the whole-tree models (462–651 kg m⁻³) was not as large as the SilviScan ADD data (320–912 kg m⁻³). For coarseness, the NIR-SilviScan calibration also had a high R² (0.86), while the NIR-Mean and NIR-Weighted models had only one factor recommended and were much weaker (R² for both=0.42). The narrow range for the FQA whole-tree values for coarseness (443–599 μg m⁻¹) compared to the SilviScan values (359–800 μg m⁻¹) may explain the inferior models. The reduced range for the whole-tree values compared to the SilviScan values arises as SilviScan values are measured at multiple locations within the tree as opposed to a single composite value.

ADD maps – 13-year-old trees

Maps showing ADD variation within 13-year-old loblolly pines based on predictions provided by the three different models are shown in Figure 2. In general, all maps are consistent with what is expected for ADD variation within loblolly pine trees (Burdon et al. 2004; Jordan et al. 2008; Schimleck et al. 2018a). They also show very similar patterns of variation; however, the NIR-Weighted spatial-interpolated map had lower estimates of ADD which gave a larger region of low ADD wood near the pith and shifted the high ADD region at the base closer to the periphery of the map. ADD probability (density) plots (Figure 3) are consistent with this observation. Figure 3 also shows that the shape of the probability plots based on predicted ADDs from the NIR-Mean and NIR-Weighted models were more similar than the plot based on the NIR-SilviScan spatial-specific model, though for ADD predictions greater than 575 kg m⁻³, the two model types compared very well.

Figure 2:

Maps showing within-tree variation of air-dry density (ADD) for loblolly pine trees aged 13 years.

Maps represent the average of 18 trees and were developed using Akima’s interpolation method and data provided by (a) spatial-specific ADD model (NIR-SilviScan), (b) interpolation-specific ADD model (NIR-Mean) and (c) interpolation-specific ADD model (NIR-Weighted).

Figure 3:

Probability (density) plots showing the variation of predicted density for trees aged 13 years.

ADD maps – 22-year-old trees

ADD maps for the 22-year-old trees (Figure 4) showed an increase in ADD at all heights from pith to bark, an observation consistent with maturation (increases in percent latewood and tracheid wall thickness) as loblolly pine trees grow older (Burdon et al. 2004). The lower predictions of ADD for the NIR-Weighted model were again apparent, while the higher ADD predictions of the NIR-Mean model were more pronounced compared to the maps for the 13-year-old trees, with a zone of very high density at the base of the tree. The probability plots (Figure 5) demonstrated the age-related shifts in ADD (also observed by Schimleck et al. 2009) and the tendency toward higher ADD predictions by the NIR-Mean model. While the maps were similar, the two spatial-interpolated models (NIR-Mean, NIR-Weighted) did not match the range of predictions provided by the spatial-specific (NIR-SilviScan) model, with the NIR-Mean and NIR-Weighted models tending to over- and under-predict the ADD, respectively.

Figure 4:

Maps showing within-tree variation of density (ADD) for loblolly pine trees aged 22 years.

Maps represent the average of 18 trees and were developed using Akima’s interpolation and data provided by (a) spatial-specific ADD model (NIR-SilviScan), (b) interpolation-specific ADD model (NIR-Mean) and (c) interpolation-specific ADD model (NIR-Weighted).

Figure 5:

Probability (density) plots showing the variation of predicted density for trees aged 22 years.

Coarseness maps – 13-year-old trees

Overall, the pattern of coarseness variation (Figure 6) was consistent with the maps for ADD, i.e. both properties demonstrate a radial increase at all heights (Burdon et al. 2004). However, the coarseness maps were not as consistent between the models as the ADD maps owing to the variable changes that occur in tracheid coarseness during pulping (Scallan and Green 1975; Kibblewhite and Evans 2001; Çöpür et al. 2005). The NIR-SilviScan coarseness map shows greater within-tree variation compared to the NIR-Mean and NIR-Weighted spatial-interpolated maps which do not have coarseness values below 425 μg m⁻¹. Probability plots (Figure 7) show NIR-SilviScan predictions had a much larger range than the NIR-Mean and NIR-Weighted approaches and confirm that the three sets of predictions were less consistent than for ADD.

Figure 6:

Maps showing within-tree variation of tracheid coarseness for loblolly pine trees aged 13 years.

Maps represent the average of 18 trees and were developed using Akima’s interpolation and data provided by (a) spatial-specific coarseness model (NIR-SilviScan), (b) interpolation-specific coarseness model (NIR-Mean) and (c) interpolation-specific coarseness model (NIR-Weighted).

Figure 7:

Probability (density) plots showing the variation of predicted coarseness for trees aged 13 years.

Coarseness maps – 22-year-old trees

Coarseness maps for the 22-year-old trees showed again the radial increase of coarseness as cambial age increases (Figure 8). Compared to the 13-year-old tree maps, there was an increase in the coarseness values at or greater than 600 μg m⁻¹ with the NIR-SilviScan map (spatial-specific) having the highest coarseness values. The probability plots (Figure 9) demonstrated the age-related shifts in coarseness. The distributions were more similar but as with the 13-year-old trees, the NIR-SilviScan model provided data with a wider range.

Figure 8:

Maps showing within-tree variation of tracheid coarseness for loblolly pine trees aged 22 years.

Maps represent the average of 18 trees and were developed using Akima’s interpolation and data provided by (a) spatial-specific ADD model (NIR-SilviScan), (b) interpolation-specific ADD model (NIR-Mean) and (c) interpolation-specific ADD model (NIR-Weighted).

Figure 9:

Probability (density) plots showing the variation of predicted coarseness for trees aged 22 years.

When the predicted values for the NIR-SilviScan model were plotted against the NIR-Mean and NIR-Weighted models (Figure 10), the datasets agreed very well for the ADD, with both models having a coefficient of determination (R²) equal to or greater than 0.90. The tendency of the NIR-Mean model to over-predict the ADD was noticeable over the range of ADD values, while the under-predictions of the NIR-Weighted model were only apparent for densities greater than 550–600 kg m⁻³. Overall differences in predictions were slightly smaller for the NIR-Weighted model [root mean square error (RMSE)=32.3 kg m⁻³ compared to 32.6 kg m⁻³]. There was a much higher deviation between the NIR-SilviScan model and the NIR-Mean and NIR-Weighted models for predicted coarseness compared to predicted density as expected, given the adjustment done on the data in an attempt to correct for differences in coarseness that occur during pulping. Both the NIR-Mean and NIR-Weighted models had the same prediction statistics (R²=0.83, 41.2 μg m⁻¹) as the NIR-SilviScan coarseness.

Figure 10:

Plots of NIR-SilviScan-predicted density vs. NIR-Mean- and NIR-Weighted- predicted density and NIR-SilviScan-predicted coarseness vs. NIR-Mean- and NIR-Weighted-predicted coarseness.

The blue line represents the line of equivalence.

Table 2 provides a comparison of measured and predicted data at the whole-tree level. For both the 13- and 22-year-old trees, the NIR-Weighted model (524 and 584 kg m⁻³) yielded a closer prediction to the measured ADD (542 and 591 kg m⁻³) than the NIR-Mean model (576 and 637 kg m⁻³). Differences between measured whole-tree ADD and NIR-SilviScan ADD (570 and 605 kg m⁻³) were also observed and can be expected because the measurements are based on different principles, one gravimetric on whole-tree chips, the other X-ray densitometry on radial strips. For coarseness, again the NIR-Weighted model (495 and 535 μg m⁻¹) was more similar to measured data (500 and 545 μg m⁻¹) than the NIR-Mean model (524 and 566 μg m⁻¹). As expected, greater differences were observed for whole-tree vs. NIR-SilviScan coarseness (488 and 588 μg m⁻¹) than what was observed for density, due to the variable impact pulping has on tracheid dimensions, and our use of the coarseness adjustment model developed for radiata pine. For the four sets of comparisons, the NIR-Weighted model was more accurate than the NIR-Mean model at estimating measured whole-tree values because weighted NIR spectra better represent the basal area of each NIR radial measurement within a disc.

Maps showing the within-tree variation of the ADD and coarseness for loblolly pine trees aged 13 and 22 years were developed using data provided by two different types of NIR-based models, referred to in this study as spatial-specific (NIR-SilviScan) and spatial-interpolated (NIR-Mean and NIR-Weighted). The two approaches fundamentally differ with the wood property data utilized measured at very different scales. For the spatial-specific model, the NIR and SilviScan data were measured in 10-mm increments, while for the spatial-interpolated models (one utilizing basal-area weighted average whole-tree spectra, and the other arithmetic average whole-tree spectra) the ADD and coarseness were measured on whole-tree composite samples.

ADD maps at both ages, obtained using the different NIR-based models, were compared and patterns of within-tree variation were similar regardless of the models used. Minor differences were observed in the maps and were confirmed by probability plots; however, the magnitude of the differences were relatively small (Figures 3 and 5). Tracheid coarseness maps (Figures 6 and 8) were also very similar and for both properties and differences are not sufficient to change the interpretation of how ADD and tracheid coarseness vary within loblolly pine trees.

The result is highly significant for future studies of wood property variation within trees. It demonstrates that properties typically measured on whole-tree composite samples (for example, pulp yield, but in this study, ADD and tracheid coarseness) can be used to develop NIR-based whole-tree models that can be used to estimate properties utilizing spectra collected at much higher spatial resolution (in this study from radial strips in 10-mm increments). Our findings also provide confidence that plots of within-sample variability based on NIR-HSI data (Thumm et al. 2010; Schimleck et al. 2018b), where whole-sample wood property data and average spectra are used to predict property variation at high pixel-by-pixel resolution, give realistic representations of within-sample variability.

The comparative approach of whole-tree maps utilized in this study was only possible as NIR spectra were measured in 10-mm increments at multiple heights and radial positions (allowing whole-tree average spectra to be obtained for each tree), and we had wood properties measured at both a specific level (radial positions, 10-mm increments) and the whole-tree level (on composite chips). If our approach were to be utilized in future studies of within-tree variation, NIR spectra would still need to be measured at multiple positions on radial samples and at multiple heights. For this study, preparation of samples and collection of NIR spectra (2569 in total) was a considerable undertaking, and now the equivalent amount of data can be quickly collected from radial strips or cross-sectional discs using NIR-HSI.

This study compared tree maps of density adjusted to air-dry conditions measured from radial strips using SilviScan and X-ray densitometry, to density adjusted to air-dry conditions measured using water displacement from chips that were processed from whole trees. Also explored were maps of coarseness measured on dry samples using SilviScan, and coarseness measured on macerated pulp. The technique explored here illustrates that the spatial-interpolated method was successful for generating information on within-tree variation and the results can be applied to other properties that are commonly determined on whole-tree chips, including wood chemistry parameters such as pulp yield, as well as handsheet properties, provided a suitable NIR model can be calculated. The spatial-interpolated approach is not suitable for some properties such as MFA, due to the difficulty in obtaining these measurements on macerated tracheids as compared to the SilviScan system using X-ray diffraction.