# Developmental Trends of Black Spruce Fibre Attributes in Maturing Plantations

Authors: Peter F. Newton
Journal: International Journal of Forestry Research
 (2016)

---

## Abstract

This study assessed the temporal developmental patterns of commercially relevant fibre attributes (tracheid length and diameters, wall thickness, specific surface area, wood density, microfibril angle, fibre coarseness, and modulus of elasticity) and their interrelationships within maturing black spruce (Picea mariana (Mill.) B.S.P.) plantations. A size-based stratified random sample procedure within 5 semimature plantations located in the Canadian Boreal Forest Region was used to select 50 trees from which radial cross-sectional xylem sequences at breast-height (1.3 m) were cut and analyzed. Statistically, the graphical and linear correlation analyses indicated that the attributes exhibited significant (p≤0.05) relationships among themselves and with morphological tree characteristics. Relative variation of each annually measured attribute declined with increasing size class (basal area quintile). The transitional shifts in temporal correlation patterns occurring at the time of approximate crown closure where suggestive of intrinsic differences in juvenile and mature wood formation processes. The temporal cumulative development patterns of all 8 of the annually measured attributes varied systematically with tree size and exhibited the most rapid rates of change before the trees reached a cambial age of 20 years. At approximately 50 years after establishment, plantation mean attribute values were not dissimilar from those reported for more mature natural-origin stands.


---

## Body

## 1. Introduction

Black spruce (Picea mariana (Mill.) B.S.P.) is a preferred reforestation species used in the establishment of plantations throughout the Canadian Boreal Forest Region [1]. Its desirable fibre characteristics and its ability to grow on a wide range of sites have combined to make black spruce an important industrial feedstock for the production of a broad array of end-products (e.g., dimensional lumber, pulp and paper, and composite wood products [2]). Black spruce plantations are an important component of the current and future industrial wood supply and are expected to provide a broad array of ecosystem services over their rotations (sensu [3]). The expectation is that these plantations will be well managed in accordance with an intensive silvicultural regime and hence will receive the required silvicultural treatments that promote end-product quality over their rotations [4]. For example, when plantations become overstocked and experience the resultant consequences of lower growth and increased mortality, treatments such as commercial thinning (CT) are expected to be carried out.Operationally, a number of underlying factors that include the establishment of narrowly spaced plantations and high levels of natural ingress have resulted in an increased occurrence of overstocked plantations across the boreal landscape. Summary statistics based on 285 permanent growth and yield sample plots established in upland black spruce plantations across boreal Ontario revealed a relatively high level of density-stress within maturing plantations (i.e., mean absolute density of 2991 stems/ha with an associated relative density index of 0.5 at an average age of 28 years [6]). Given that a relative density index of 0.5 is the approximate threshold value at which imminent competition mortality or self-thinning is likely to occur [6], many plantations are likely to experience the adverse effects of increased competition which could include the loss of merchantable-sized trees. In response, forest managers are actively revising their initial crop plans in terms of adding density control treatments, such as CT [7]. Additionally, apart from ameliorating the effects of density-stress, thinning yields (residues) arising from CT treatments potentially represent an important interim supplemental source of fibre at midrotation (e.g., [8]).Among the principal determinates governing the end-product potential of harvested trees are their internal fibre characteristics. For example, fibre attributes, such as tracheid dimensions (length and radial and tangential diameters (fl, dr, and dt, resp.)), tracheid wall thickness (wt), specific surface area (sa), density (wd), microfibril angle (ma or MFA), fibre coarseness (co), and modulus of elasticity (me or MOE), directly influence end-product quality and quantity (Table 1 (sensu [5])). Consequently, attaining an understanding of the status of these commercially relevant fibre attributes as plantations develop and identifying factors that may influence their development could provide the prerequisite knowledge for forecasting end-product potential, in addition to helping guide forest managers in designing optimal silvicultural treatment schedules.Table 1
Product-based performance measures and their relationship with fibre attributes (sensu [5]).
Product category Performance measure Functional relationship with fibre attribute Pulp and paper Tensile strength ∝ fibre length, (wall thickness)−1 Tear strength ∝ fibre length, coarseness Stretch ∝ microfibril angle Bulk ∝ wall thickness, (fibre width)−1 Light scattering ∝ (wall thickness)−1 Sheet formation ∝ (fibre length)−1 Collapsibility ∝ wall thickness Yield ∝ density Solid wood and composites Strength ∝ density, (microfibril angle)−1 Stiffness ∝ density, modulus of elasticity, (microfibril angle)−1Thus, the objectives of this study were to investigate static and dynamic correlative patterns of association and the temporal developmental trends of commercially relevant fibre attributes within maturing black spruce plantations. The overall goals of this study were to provide inferences on the nature of the formulation pattern of these attributes, determine if developmental trends varied by tree size (i.e., hierarchical position within the stand as measured by basal area quintile), and infer potential end-products of conceptual CT yields based on the observed attribute values at midrotation relative to those observed within older natural-origin stands.


## 2. Materials and Methods

### 2.1. Data Acquisition and Processing

Breast-height (1.3 m) cross-sectional disks were obtained from 50 plantation-grown sample trees during the autumn of 2006 and 2007. These trees were a sample subset taken from a larger investigation studying tree taper and diameter growth patterns of boreal conifers (e.g., see [9–11] for specifics). The 50 trees were sampled from 5 even-aged monospecific black spruce plantations that were situated on upland sites throughout Forest Sections B9 (Superior) and B11 (Upper English River) of the Canadian Boreal Forest Region [1]. Geographically, the plantations were located within a rectangular region bounded by latitudes of 48° N on the south and 50° N on the north and longitudes of 88° W on the east and 92° W on the west (i.e., specific latitude/longitudes of the plantations were 49°05′12′′/88°09′58′′ (Plantation code W15); 48°56′06′′/90°35′51′′ (W18); 49°57′52′′/92°29′39′′ (W19); 50°11′′47′/91°42′16′′ (W20); and 48°42′11′′/90°08′06′′ (W23)). All five plantations shared a common silvicultural history: planted at approximately 2500 stems/ha (2 × 2 m intertree spacing) in the late 1960s or early 1970s on sites that were mechanically scarified.The plantations varied slightly in terms of their productivity (medium-low to medium-high site qualities as inferred from their site index values (mean/min/max values of 16.4/14.8/18.3 m at 50 years, resp. [12]), ages (mean/min/max breast-height ages of 36/34/43, resp.), absolute stand densities (mean/min/max densities (stems/ha) of 3793/2076/5579, resp.) and site occupancy (fully stocked to overstocked as measured by relative density index (mean/min/max values of 0.76/0.51/1.02, resp. [13]). Variation in densities among the plantations was largely attributed to differential rates in natural ingress which occurred immediately following establishment.Three variable-sized circular temporary sample plots were established within each plantation. The initial size of the sample plots was set at 400 m2 which was, if necessary, increased by 100 m2 increments until a minimum of 80 trees per plot were included. All living black spruce trees were sequentially numbered and measured for diameter at breast-height (1.3 m) outside-bark (D; cm) within each plot. The trees were then stratified into basal area quintiles from which 1 tree, which was classified as planted and did not exhibit any visible deformities, such as forks, major stem injuries, and dead or broken tops, was selected for destructive stem analysis. This stratified random sampling protocol ensured that all size and crown classes were represented and indirectly provided an analytical framework for assessing fibre attribute variation and developmental differences throughout the stand structure or size hierarchy. Subsequent to felling each sample tree, total height (H; m), breast-height age (A; yr), live crown ratio (Cr (%), calculated as the length of the live crown divided by total tree height and expressed as a percentage), and diameters of the largest dead (cm) and living (cm) branches were measured. Table 2 provides a summary of the principal mensurational characteristics of the sampled plantations and the individual sample trees selected.Table 2
Descriptive statistical summary of the mensurational characteristics of the 5 plantations and 47 sample trees.
Variable Symbol Unit Mean Standard error Minimum Maximum Plantation level Age at breast-height of dominant trees  yr 36 2 34 43 Quadratic mean diameter  cm 11.0 1.45 6.4 15.8 Mean dominant height  m 12.9 0.90 10.5 16.5 Site index  m 16.4 0.70 14.8 18.3 Density  stems/ha 3793 449 2076 5579 Basal area  m2/ha 30.0 4.78 17.8 45.4 Total volume  m3/ha 189.1 43.27 87.1 348.3 Merchantable volume  m3/ha 125.9 44.97 18.6 292.8 Mean volume per tree  dm3 60.2 21.30 15.6 146.4 Relative density index  %/100 0.76 0.075 0.51 1.02 Tree level Breast-height age A yr 33 4 25 43 Height H m 10.9 1.03 6.7 17.9 Diameter at breast-height D cm 12.2 1.43 5.4 24.8 Maximum dead branch diameter  cm 1.3 0.09 0.5 2.2 Maximum live branch diameter  cm 1.5 0.10 0.1 2.5 Live crown ratio C r % 59.4 4.91 37.2 96.7The destructive stem analysis protocol consisted of felling each tree at stump height (approximately 30 cm above the ground surface) and obtaining 13 cross-sectional samples along the stem, specifically at 0.30, 0.5, 0.9, and 1.3 m, and thereafter at every 10% height interval based on remaining distance between breast-height and the stem tip [9]. These cross-sectional samples were then placed in burlap bags, transported to the research laboratory, and stored at –10°C until 24 hours prior to processing. For the purposes of this study, 2 breast-height cross-sectional samples were randomly selected from each of the basal area quintiles within each plantation. A bark-to-pith-to-bark 2 cm square transverse sample along the geometric mean diameter was sawn from each disk, yielding a total of 50 cross-sectional transverse samples (1 transverse breast-height sample/tree × 2 trees/quintile × 5 quintiles/plantation × 5 plantations). The annual ring-width sequence and age of each sample were assessed for age inconsistencies and the presence of missing or partial rings using a binocular microscope. Age inconsistencies were found in 3 of the samples for which the observed age exceeded the age of the plantation. This finding would suggest that these trees originated from nonplanted residual saplings that were present in the previous stand. Consequently, these 3 trees were removed from the analysis. One pith-to-bark radial sequence was then sawn from each of the transverse samples. This resulted in a total of 47 annual-ring-width radial sequences available for fibre attribute determination.The anatomical characterization of the fibre attributes along each annual-ring-width sequence was carried out using the SilviScan-3 system [14] and the Fibre Quality Analyzer [15]. Briefly, SilviScan-3 is the latest iteration of the integrated wood analysis system originally developed by Dr. Robert Evans and his research colleagues at the CSIRO’s (Commonwealth Scientific and Industrial Research Organisation, Australia) Forestry and Forest Products Division. The system combines automatic image acquisition and analysis (cell scanner), X-ray densitometry, and X-ray diffractometry, to determine a multitude of fibre characteristics. For the purposes of this study, the following attributes were used from the Silviscan-3 output: (1) fibre dimensions inclusive of radial and tangential diameters, wall thickness, and specific surface area as determined via X-ray densitometry [16]; (2) fibre coarseness and density as determined via X-ray densitometry [16]; (3) microfibril angle as determined via X-ray diffraction [17]; and (4) modulus of elasticity as determined from a combination of X-ray densitometry and diffraction measurements [18]. Note that the samples were extracted with acetone in order to remove resins which may influence the density estimates. Specifically, the samples were soaked in acetone for 12 hours and then extracted for 8 hours at 70°C using a modified Soxhlet extraction system. After extraction, the samples were air-dried for approximately 12 hours and then conditioned to a 40% relative humidity at 20°C.Separately, mean fibre length was determined for the entire sequence using a high resolution Fibre Quality Analyzer (OpTest Equipment Inc.) [15, 19] employing the following methodology (Source: Dr. Tong, FPInnovations Inc., Vancouver, British Columbia, Canada). The transverse radial sequence was sectioned into 2 age class segments for samples less than 30 years of age and into 3 age class segments for samples greater than or equal to 30 years. The high resolution FQA, which is comprised of hydraulic, optical, and image processing systems, was used to derive a mean fibre length estimate from a composite (bulk) macerated sample consisting of fibres from all the segments combined. More precisely, a dilute suspension of pulp fibres from each segment was obtained and processed through optical and imaging systems by a flow cell that orients fibres hydrodynamically into a 2-dimensional plane. The digital imaging system acquired and analyzed the images of the oriented fibres from which their length and other morphological properties were determined. To obtain the pulp sample, each segment was boiled in deionized water for approximately 4 hours after which the temperature was reduced to 70°C. Subsequently, the samples were soaked for an additional 12 hours. They were then macerated in a solution of hydrogen peroxide (35%) and glacial acetic acid (1 : 1 v/v) for 48 hours at a temperature of 70°C. The obtained pulp was washed, dispersed using a mixer, and filtered using a 150 mesh screen. Fibre length was then determined by sampling and weighing a predetermined amount of the dry pulp from the resulting handsheet, diluting it to a predetermined consistency, and analyzing the sample (in duplicate) using the FQA. Based on the number of fibres in each length class and mean fibre length within each class, an overall mean weight-weighted fibre length was calculated. Table 3 provides a statistical summary of the fibre attributes of the selected sample trees.Table 3
Statistical summary of the area-weighted fibre attribute values derived from the 47 breast-height cross-sectional sample disks at the time of sampling.
Variable Symbol (nominal equivalent) Unit Mean Standard error Minimum Maximum CVa (%) Fibre length (weight-weighted) f l mm 2.4 0.2 1.8 2.8 8.8 Wood density w d kg/m3 487.3 42.5 399.7 574.3 8.7 Microfibril angle m a (MFA) ° 15.0 3.6 9.4 23.6 23.9 Modulus of elasticity m e (MOE) GPa 13.0 2.7 8.5 18.0 20.4 Coarseness c o µg/m 332.4 31.2 249.7 379.4 9.4 Wall thickness w t µm 2.3 0.2 1.9 2.6 7.9 Tracheid radial diameter d r µm 27.6 2.0 22.7 33.5 7.1 Tracheid tangential diameter d t µm 25.9 1.4 22.6 29.0 5.4 Specific surface area s a m2/kg 336.1 25.0 305.0 400.8 7.4 aCoefficient of variation.


### 2.2. Computations and Analyses

#### 2.2.1. Correlative Relationships

Graphical analyses were used to determine the most appropriate linear expression for investigating (1) bivariate correlative relationships among attributes at the time of sampling (i.e., pairwise linear association among fl, wd, ma, me, co, wt, dr, dt, and sa); (2) bivariate correlative relationships between attributes (fl, wd, ma, me, co, wt, dr, dt, and sa) at the time of sampling and individual-tree morphological metrics (D,H,A, and Cr); and (3) temporal annual patterns of bivariate correlation between attributes (i.e., pairwise linear association among ring width area (ar (mm2)), wd, ma, me, co, wt, dr, dt, and sa by cambial age). Pairwise linearity was verified graphically using bivariate scatterplots and then statistically quantified via the Pearson moment correlation coefficient (r).Attribute variation among the sample trees at the time of sampling was measured using the coefficient of variation. Specifically, the ring area-weighted mean and associated standard derivation were calculated forwd, ma, me, co, wt, dr, dt, and sa along each pith-to-bark radial sequence, from which the coefficient of variation was computed: (1)cvk=skmk=sk2∑j=1Jaj2/∑j=1Jaj21/2∑j=1Jvkjaj/∑j=1Jaj,where cvk, sk, sk2, and mk are, respectively, the area-weighted sequence-specific coefficient of variation, standard deviation, variance, and mean value of thekth fibre attribute, aj is the area of the annual ring (mm2) corresponding to thejth cambial age (j=1,…,J), and vkj is the mean annual ring-width-area weighed value specific to thekth attribute andjth cambial age.


#### 2.2.2. Accumulative Computations and Resultant Indices

For each sequence, the cumulative annual-ring-area-weighted moving average was calculated for each attribute (wd, ma, me, co, wt, dr, dt, and sa) in the pith-to-bark direction: (2)Vki=∑j=1ivkjaj∑j=1iaj,where Vki is the cumulative area-weighted moving average value specific to thekth attribute calculated to the ith cambial age. Note that Wd, Ma, Me, Co, Wt, Dr, Dt, and Sa are used to denote the cumulative moving average value (i.e., Vki in (2)) for attributes wd, ma, me, co, wt, dr, dt, and sa, respectively. The cumulative moving average value reflects the accumulated status of a given attribute up to the ith cambial age and hence is indicative of the overall product potential of a tree if harvested at that specific age. Three-dimensional graphics were utilized in evaluating and summarizing the pattern of cumulative development of each variable over time by tree size class (basal area quintile). Note that Fortran-based analytical programs were written in order to calculate the cumulative values, evaluate and select the most appropriate bivariate correlative relationships, and generate the descriptive statistics. The graphical analysis and the statistical computations and testing were carried out employing Statistica (V12; Dell Inc.).


## 2.1. Data Acquisition and Processing

Breast-height (1.3 m) cross-sectional disks were obtained from 50 plantation-grown sample trees during the autumn of 2006 and 2007. These trees were a sample subset taken from a larger investigation studying tree taper and diameter growth patterns of boreal conifers (e.g., see [9–11] for specifics). The 50 trees were sampled from 5 even-aged monospecific black spruce plantations that were situated on upland sites throughout Forest Sections B9 (Superior) and B11 (Upper English River) of the Canadian Boreal Forest Region [1]. Geographically, the plantations were located within a rectangular region bounded by latitudes of 48° N on the south and 50° N on the north and longitudes of 88° W on the east and 92° W on the west (i.e., specific latitude/longitudes of the plantations were 49°05′12′′/88°09′58′′ (Plantation code W15); 48°56′06′′/90°35′51′′ (W18); 49°57′52′′/92°29′39′′ (W19); 50°11′′47′/91°42′16′′ (W20); and 48°42′11′′/90°08′06′′ (W23)). All five plantations shared a common silvicultural history: planted at approximately 2500 stems/ha (2 × 2 m intertree spacing) in the late 1960s or early 1970s on sites that were mechanically scarified.The plantations varied slightly in terms of their productivity (medium-low to medium-high site qualities as inferred from their site index values (mean/min/max values of 16.4/14.8/18.3 m at 50 years, resp. [12]), ages (mean/min/max breast-height ages of 36/34/43, resp.), absolute stand densities (mean/min/max densities (stems/ha) of 3793/2076/5579, resp.) and site occupancy (fully stocked to overstocked as measured by relative density index (mean/min/max values of 0.76/0.51/1.02, resp. [13]). Variation in densities among the plantations was largely attributed to differential rates in natural ingress which occurred immediately following establishment.Three variable-sized circular temporary sample plots were established within each plantation. The initial size of the sample plots was set at 400 m2 which was, if necessary, increased by 100 m2 increments until a minimum of 80 trees per plot were included. All living black spruce trees were sequentially numbered and measured for diameter at breast-height (1.3 m) outside-bark (D; cm) within each plot. The trees were then stratified into basal area quintiles from which 1 tree, which was classified as planted and did not exhibit any visible deformities, such as forks, major stem injuries, and dead or broken tops, was selected for destructive stem analysis. This stratified random sampling protocol ensured that all size and crown classes were represented and indirectly provided an analytical framework for assessing fibre attribute variation and developmental differences throughout the stand structure or size hierarchy. Subsequent to felling each sample tree, total height (H; m), breast-height age (A; yr), live crown ratio (Cr (%), calculated as the length of the live crown divided by total tree height and expressed as a percentage), and diameters of the largest dead (cm) and living (cm) branches were measured. Table 2 provides a summary of the principal mensurational characteristics of the sampled plantations and the individual sample trees selected.Table 2
Descriptive statistical summary of the mensurational characteristics of the 5 plantations and 47 sample trees.
Variable Symbol Unit Mean Standard error Minimum Maximum Plantation level Age at breast-height of dominant trees  yr 36 2 34 43 Quadratic mean diameter  cm 11.0 1.45 6.4 15.8 Mean dominant height  m 12.9 0.90 10.5 16.5 Site index  m 16.4 0.70 14.8 18.3 Density  stems/ha 3793 449 2076 5579 Basal area  m2/ha 30.0 4.78 17.8 45.4 Total volume  m3/ha 189.1 43.27 87.1 348.3 Merchantable volume  m3/ha 125.9 44.97 18.6 292.8 Mean volume per tree  dm3 60.2 21.30 15.6 146.4 Relative density index  %/100 0.76 0.075 0.51 1.02 Tree level Breast-height age A yr 33 4 25 43 Height H m 10.9 1.03 6.7 17.9 Diameter at breast-height D cm 12.2 1.43 5.4 24.8 Maximum dead branch diameter  cm 1.3 0.09 0.5 2.2 Maximum live branch diameter  cm 1.5 0.10 0.1 2.5 Live crown ratio C r % 59.4 4.91 37.2 96.7The destructive stem analysis protocol consisted of felling each tree at stump height (approximately 30 cm above the ground surface) and obtaining 13 cross-sectional samples along the stem, specifically at 0.30, 0.5, 0.9, and 1.3 m, and thereafter at every 10% height interval based on remaining distance between breast-height and the stem tip [9]. These cross-sectional samples were then placed in burlap bags, transported to the research laboratory, and stored at –10°C until 24 hours prior to processing. For the purposes of this study, 2 breast-height cross-sectional samples were randomly selected from each of the basal area quintiles within each plantation. A bark-to-pith-to-bark 2 cm square transverse sample along the geometric mean diameter was sawn from each disk, yielding a total of 50 cross-sectional transverse samples (1 transverse breast-height sample/tree × 2 trees/quintile × 5 quintiles/plantation × 5 plantations). The annual ring-width sequence and age of each sample were assessed for age inconsistencies and the presence of missing or partial rings using a binocular microscope. Age inconsistencies were found in 3 of the samples for which the observed age exceeded the age of the plantation. This finding would suggest that these trees originated from nonplanted residual saplings that were present in the previous stand. Consequently, these 3 trees were removed from the analysis. One pith-to-bark radial sequence was then sawn from each of the transverse samples. This resulted in a total of 47 annual-ring-width radial sequences available for fibre attribute determination.The anatomical characterization of the fibre attributes along each annual-ring-width sequence was carried out using the SilviScan-3 system [14] and the Fibre Quality Analyzer [15]. Briefly, SilviScan-3 is the latest iteration of the integrated wood analysis system originally developed by Dr. Robert Evans and his research colleagues at the CSIRO’s (Commonwealth Scientific and Industrial Research Organisation, Australia) Forestry and Forest Products Division. The system combines automatic image acquisition and analysis (cell scanner), X-ray densitometry, and X-ray diffractometry, to determine a multitude of fibre characteristics. For the purposes of this study, the following attributes were used from the Silviscan-3 output: (1) fibre dimensions inclusive of radial and tangential diameters, wall thickness, and specific surface area as determined via X-ray densitometry [16]; (2) fibre coarseness and density as determined via X-ray densitometry [16]; (3) microfibril angle as determined via X-ray diffraction [17]; and (4) modulus of elasticity as determined from a combination of X-ray densitometry and diffraction measurements [18]. Note that the samples were extracted with acetone in order to remove resins which may influence the density estimates. Specifically, the samples were soaked in acetone for 12 hours and then extracted for 8 hours at 70°C using a modified Soxhlet extraction system. After extraction, the samples were air-dried for approximately 12 hours and then conditioned to a 40% relative humidity at 20°C.Separately, mean fibre length was determined for the entire sequence using a high resolution Fibre Quality Analyzer (OpTest Equipment Inc.) [15, 19] employing the following methodology (Source: Dr. Tong, FPInnovations Inc., Vancouver, British Columbia, Canada). The transverse radial sequence was sectioned into 2 age class segments for samples less than 30 years of age and into 3 age class segments for samples greater than or equal to 30 years. The high resolution FQA, which is comprised of hydraulic, optical, and image processing systems, was used to derive a mean fibre length estimate from a composite (bulk) macerated sample consisting of fibres from all the segments combined. More precisely, a dilute suspension of pulp fibres from each segment was obtained and processed through optical and imaging systems by a flow cell that orients fibres hydrodynamically into a 2-dimensional plane. The digital imaging system acquired and analyzed the images of the oriented fibres from which their length and other morphological properties were determined. To obtain the pulp sample, each segment was boiled in deionized water for approximately 4 hours after which the temperature was reduced to 70°C. Subsequently, the samples were soaked for an additional 12 hours. They were then macerated in a solution of hydrogen peroxide (35%) and glacial acetic acid (1 : 1 v/v) for 48 hours at a temperature of 70°C. The obtained pulp was washed, dispersed using a mixer, and filtered using a 150 mesh screen. Fibre length was then determined by sampling and weighing a predetermined amount of the dry pulp from the resulting handsheet, diluting it to a predetermined consistency, and analyzing the sample (in duplicate) using the FQA. Based on the number of fibres in each length class and mean fibre length within each class, an overall mean weight-weighted fibre length was calculated. Table 3 provides a statistical summary of the fibre attributes of the selected sample trees.Table 3
Statistical summary of the area-weighted fibre attribute values derived from the 47 breast-height cross-sectional sample disks at the time of sampling.
Variable Symbol (nominal equivalent) Unit Mean Standard error Minimum Maximum CVa (%) Fibre length (weight-weighted) f l mm 2.4 0.2 1.8 2.8 8.8 Wood density w d kg/m3 487.3 42.5 399.7 574.3 8.7 Microfibril angle m a (MFA) ° 15.0 3.6 9.4 23.6 23.9 Modulus of elasticity m e (MOE) GPa 13.0 2.7 8.5 18.0 20.4 Coarseness c o µg/m 332.4 31.2 249.7 379.4 9.4 Wall thickness w t µm 2.3 0.2 1.9 2.6 7.9 Tracheid radial diameter d r µm 27.6 2.0 22.7 33.5 7.1 Tracheid tangential diameter d t µm 25.9 1.4 22.6 29.0 5.4 Specific surface area s a m2/kg 336.1 25.0 305.0 400.8 7.4 aCoefficient of variation.


## 2.2. Computations and Analyses

### 2.2.1. Correlative Relationships

Graphical analyses were used to determine the most appropriate linear expression for investigating (1) bivariate correlative relationships among attributes at the time of sampling (i.e., pairwise linear association among fl, wd, ma, me, co, wt, dr, dt, and sa); (2) bivariate correlative relationships between attributes (fl, wd, ma, me, co, wt, dr, dt, and sa) at the time of sampling and individual-tree morphological metrics (D,H,A, and Cr); and (3) temporal annual patterns of bivariate correlation between attributes (i.e., pairwise linear association among ring width area (ar (mm2)), wd, ma, me, co, wt, dr, dt, and sa by cambial age). Pairwise linearity was verified graphically using bivariate scatterplots and then statistically quantified via the Pearson moment correlation coefficient (r).Attribute variation among the sample trees at the time of sampling was measured using the coefficient of variation. Specifically, the ring area-weighted mean and associated standard derivation were calculated forwd, ma, me, co, wt, dr, dt, and sa along each pith-to-bark radial sequence, from which the coefficient of variation was computed: (1)cvk=skmk=sk2∑j=1Jaj2/∑j=1Jaj21/2∑j=1Jvkjaj/∑j=1Jaj,where cvk, sk, sk2, and mk are, respectively, the area-weighted sequence-specific coefficient of variation, standard deviation, variance, and mean value of thekth fibre attribute, aj is the area of the annual ring (mm2) corresponding to thejth cambial age (j=1,…,J), and vkj is the mean annual ring-width-area weighed value specific to thekth attribute andjth cambial age.


### 2.2.2. Accumulative Computations and Resultant Indices

For each sequence, the cumulative annual-ring-area-weighted moving average was calculated for each attribute (wd, ma, me, co, wt, dr, dt, and sa) in the pith-to-bark direction: (2)Vki=∑j=1ivkjaj∑j=1iaj,where Vki is the cumulative area-weighted moving average value specific to thekth attribute calculated to the ith cambial age. Note that Wd, Ma, Me, Co, Wt, Dr, Dt, and Sa are used to denote the cumulative moving average value (i.e., Vki in (2)) for attributes wd, ma, me, co, wt, dr, dt, and sa, respectively. The cumulative moving average value reflects the accumulated status of a given attribute up to the ith cambial age and hence is indicative of the overall product potential of a tree if harvested at that specific age. Three-dimensional graphics were utilized in evaluating and summarizing the pattern of cumulative development of each variable over time by tree size class (basal area quintile). Note that Fortran-based analytical programs were written in order to calculate the cumulative values, evaluate and select the most appropriate bivariate correlative relationships, and generate the descriptive statistics. The graphical analysis and the statistical computations and testing were carried out employing Statistica (V12; Dell Inc.).


## 2.2.1. Correlative Relationships

Graphical analyses were used to determine the most appropriate linear expression for investigating (1) bivariate correlative relationships among attributes at the time of sampling (i.e., pairwise linear association among fl, wd, ma, me, co, wt, dr, dt, and sa); (2) bivariate correlative relationships between attributes (fl, wd, ma, me, co, wt, dr, dt, and sa) at the time of sampling and individual-tree morphological metrics (D,H,A, and Cr); and (3) temporal annual patterns of bivariate correlation between attributes (i.e., pairwise linear association among ring width area (ar (mm2)), wd, ma, me, co, wt, dr, dt, and sa by cambial age). Pairwise linearity was verified graphically using bivariate scatterplots and then statistically quantified via the Pearson moment correlation coefficient (r).Attribute variation among the sample trees at the time of sampling was measured using the coefficient of variation. Specifically, the ring area-weighted mean and associated standard derivation were calculated forwd, ma, me, co, wt, dr, dt, and sa along each pith-to-bark radial sequence, from which the coefficient of variation was computed: (1)cvk=skmk=sk2∑j=1Jaj2/∑j=1Jaj21/2∑j=1Jvkjaj/∑j=1Jaj,where cvk, sk, sk2, and mk are, respectively, the area-weighted sequence-specific coefficient of variation, standard deviation, variance, and mean value of thekth fibre attribute, aj is the area of the annual ring (mm2) corresponding to thejth cambial age (j=1,…,J), and vkj is the mean annual ring-width-area weighed value specific to thekth attribute andjth cambial age.


## 2.2.2. Accumulative Computations and Resultant Indices

For each sequence, the cumulative annual-ring-area-weighted moving average was calculated for each attribute (wd, ma, me, co, wt, dr, dt, and sa) in the pith-to-bark direction: (2)Vki=∑j=1ivkjaj∑j=1iaj,where Vki is the cumulative area-weighted moving average value specific to thekth attribute calculated to the ith cambial age. Note that Wd, Ma, Me, Co, Wt, Dr, Dt, and Sa are used to denote the cumulative moving average value (i.e., Vki in (2)) for attributes wd, ma, me, co, wt, dr, dt, and sa, respectively. The cumulative moving average value reflects the accumulated status of a given attribute up to the ith cambial age and hence is indicative of the overall product potential of a tree if harvested at that specific age. Three-dimensional graphics were utilized in evaluating and summarizing the pattern of cumulative development of each variable over time by tree size class (basal area quintile). Note that Fortran-based analytical programs were written in order to calculate the cumulative values, evaluate and select the most appropriate bivariate correlative relationships, and generate the descriptive statistics. The graphical analysis and the statistical computations and testing were carried out employing Statistica (V12; Dell Inc.).


## 3. Results and Discussion

### 3.1. Static Correlative Relationships

The correlation patterns between the fibre attributes largely followed expectation (Table4). Increasing fibre lengths were associated with enhanced lumber related attributes such as lower ma and increased me (stiffness). Decreases in tracheid diameters were associated with increasing wood density and stiffness as evident by the inverse dr–wd-1, dr–me-1, dt–wd-1, and dt–me-1 correlative relationships. Increasing wall thickness which is associated with increased wood density and wood stiffness was also reflected in the results: that is, the positive wt–wd and me–wt correlations. Strong correlations (r≥0.70) were observed among (1) wood density and the modulus of elasticity which are both related to solid wood products and lumber grades (wd–me) and (2) fibre dimensional attributes associated with pulp yields and paper end-products (fl–co, co–wt, co–sa-1, wt–sa-1, and dr–dt). Essentially, these results suggest that increasing wood density is accompanied by a suite of other favourable attributes which could translate into greater solid wood product yields and grades (Table 1). Similarly, long fibre lengths are accompanied with a suite of attributes which collectively are associated with enhanced pulp and paper yields and quality (Table 1).Table 4
Bivariate linear association among fibre attribute values at midrotation as measured by the Pearson product moment correlation coefficient.
Attributea Attributea f l w d m a m e c o w t d r d t f l 1.0 w d 0.0400 1.0 m a −0.3343∗ −0.4840∗ 1.0 m e 0.2995∗ 0.7765∗ −0.8981∗ 1.0 c o 0.6992∗ 0.1969∗ −0.1569 0.2044 1.0 w t 0.4476∗ 0.8201∗ −0.4429∗ 0.6765∗ 0.7201∗ 1.0 d r 0.4799∗ −0.5952∗ 0.2820 −0.4477∗ 0.5913∗ −0.0769 1.0 d t 0.5542∗ −0.4677∗ 0.1585 −0.3037∗ 0.6700∗ 0.0609 0.7072∗ 1.0 s a −0.5090∗ −0.7164∗ 0.3791∗ −0.5812∗ −0.8103∗ −0.9743∗ 0.0560 −0.2047 Note:∗significant correlation at the 0.05 probability level. aDenotations are defined in Table 3.Linear bivariate graphical patterns and significant (p≤0.05) correlations existed between some of the fibre attributes and tree level variables (Table 5). Fibre length was positively correlated with individual-tree diameter, height, and age. Tracheid length is associated with the degree of fibre-bonding and hence is directly related to the tensile and tear strength of derived paper products [20]. Consequently, the observed correlations for fibre length, suggest that the suitability of these plantation black spruce trees for pulp and paper production is enhanced with increasing tree size and age. Similar patterns of correlation between fibre attributes and diameter were observed for height (fl, wd-1, ma, me-1, dr, and dt-H). Positive correlations between age and co, dr, and dt were also detected. Attributes ma and sa were negatively correlated with live crown ratio whereas me was positively correlated. Overall, the attributes were most strongly and frequently correlated the individual-tree size metrics (diameter and height).Table 5
Bivariate linear association between tree level variables and fibre attribute values at midrotation as measured by the Pearson product moment correlation coefficient.
Attributea Tree level variablesb D H A C r f l 0.3636∗ 0.4017∗ 0.4117∗ 0.1808 w d −0.6089∗ −0.4977∗ −0.2559 0.2540 m a 0.3757∗ 0.3337∗ 0.2293 −0.4099∗ m e −0.5144∗ −0.4068∗ −0.2499 0.3752∗ c o 0.2493 0.1898 0.2946∗ 0.1878 w t −0.2832 −0.2330 0.0025 0.2820 d r 0.7186∗ 0.5728∗ 0.4377∗ −0.0671 d t 0.4764∗ 0.3795∗ 0.3477∗ −0.0096 s a 0.1785 0.1704 −0.0296 −0.3185∗ Note:∗significant correlation at the 0.05 probability level. aDenotations are defined in Table 3. bDenotations are defined in Table 2.


### 3.2. Dynamic Correlative Relationships

Figure1 graphically illustrates the pairwise correlation among the attributes on an annual basis (cambial age) as the plantations developed. Six general trends emerged: (1) noncontinuous inverse correlative patterns for the ar–wd, ar–me, ar–wt, wd–dr, me–dr,me–dt,me–sa,and dt–sa relationships; (2) noncontinuous directly proportional patterns for the ar–ma, ar–sa, wd–co, ma–dr, me–wt, co–dr, and co–dt relationships; (3) approximately constant and continuous positive correlations across all cambial ages for the wd–wt and co–wt relationships; (4) approximately constant and continuous negative correlations across all cambial ages for the wd–sa, ma–me, co–sa, and wt–sa relationships; (5) shifting negative-to-positive correlations for the ar–co, ar–dr, ar–dt, wd–me, ma–sa,and dr–dt relationships; and (6) shifting positive-to-negative correlations for the wd–ma, wd–dt, and ma–wt relationships.Figure 1
Matrix of temporal annual correlation patterns for bivariate relationships between individual attributes by cambial age. Note: (1) only Pearson moment correlation coefficients which were significant (p≤0.05) are displayed, and (2) ar, wd, ma, me, co, wt, dr, dt, and sa denote annual-ring-width area, wood density, microfibril angle, modulus of elasticity, fibre coarseness, tracheid wall thickness, radial tracheid diameter, tangential tracheid diameter, and specific surface area, respectively.These results suggest the presence of a set of developmental (age independent) invariant bivariate associations which included thewd–wt, wd–sa, ma–me, co–wt, co–sa, and wt–sa relationships. The underlying morphological (functional) relationships that exist among the attributes may be partially responsible for these results (e.g., wd=fwt). Furthermore, these temporal invariant relationships provide an opportunity to employ surrogate variables for estimating the values of other attributes, through a sequential set of nested bivariate relationships: wt=fwd→sa=fwt→co=fsa.The shifting positive-to-negative and negative-to-positive correlative patterns are plausible transitions points as the plantations reached the crown closure stage of development. This is the stage of stand development that is commonly associated with the cessation of juvenile or crown-based wood formation. These transition points coincided with cambial ages ranging from approximately 8 to 13 years as inferred graphically from Figure1. Assuming a 5-year developmental period is required for the sample trees to reach breast-height (1.3 m), this cambial age range corresponds to a total stand age of approximately 13–18 years. This total age range is similar to the juvenile-mature wood transition age that has previously been reported in the literature for black spruce (i.e., 11–21 year range [21, 22]) and is also similar to that reported for white spruce (Picea glauca (Moench) Voss) [23]. Based on predictions derived from a structural stand density management diagram for these sampled plantations [6], the observed transition age range corresponds to the approximate time of crown closure.Furthermore, the temporal positive-to-negative and negative-to-positive shifts in the correlations may indicate intrinsic differences in the underlying wood formation processes during the juvenile and postjuvenile growth stages. The majority of the attributes which exhibited transitions in their correlative patterns were directly or indirectly associated with cell geometry; for example, increasing radial diameters were associated with decreasing tangential diameters during the juvenile phase whereas they were directly proportional to each other during the maturing phase. Results derived from other studies are in accordance with this inference. For example, Lenz et al. [23] suggested that the accelerated growth which is likely the outcome of frequent anticline divisions of cambial initials partially underlies the unique geometrical characteristics of cells produced during the juvenile phase.Collectively, these results provide confirmatory support for previous research in which cambial age has been found to be associated with unique fibre formation patterns [24]. For example, microfibril angle has been used as an investigative tool to differentiate between juvenile and mature wood types [25]. The higher growth rate and production of juvenile wood during the precrown closure stage may also elicit a temporary decline in wood quality in terms of end-product potential (e.g., lower wd and me and higher ma (Table 1)). Overall, the results and supporting references suggest that physiologically induced changes in allometric relationships were partially responsible for the observed shift in correlation patterns as the black spruce plantations transitioned from the precrown closure stage of development to the postcrown closure stage of development.


### 3.3. Temporal Developmental Patterns

The temporal development trends of the cumulative area-weighted moving average values by tree size (basal area quintile) for wood density, microfibril angle, modulus of elasticity, coarseness, tracheid wall thickness, tracheid radial diameter, tracheid tangential diameter, and specific surface area are graphically presented in Figure2. Irrespective of attribute, the patterns revealed that the most rapid period of change was observed within the first decade of growth and that the rate of change varied inversely with tree size (e.g., rate of change declined with increased tree size). Among the attributes, Wd, Co, Wt, Dr, Dt, and Sa exhibited the greatest rate of change during this period whereas Ma and Me exhibited the least. Following this initial period of rapid change, (1) Wd, Co, Wt, Dr, Dt, and Sa exhibited a period of relative stability that was characterized by a marginal but continuous increase (Wd, Co, Wt, Dr, and Dt) or decrease (Sa) with increasing cambial age, and (2) Ma and Me exhibited a continuous declining (Ma) or increasing (Me) trend with increasing cambial age. For a given cambial age, Wd, Me, Co, and Wt exhibited size-dependency in which values declined with increasing basal area quintile. After approximately 40 years, the rate of change within a given basal area quintile was minimal.Figure 2
Three-dimensional visualization of mean temporal developmental trends illustrating patterns of change in attribute-specific cumulative area-weighted moving average values with increasing cambial age and size quintile. Note:Wd, Ma, Me, Co, Wt, Dr, Dt, and Sa denote the mean value up to the ith cambial age for wood density, microfibril angle, modulus of elasticity, fibre coarseness, tracheid wall thickness, radial tracheid diameter, tangential tracheid diameter, and specific surface area, respectively.The results indicated that the most rapid change within the trajectories occurred before the trees reached a cambial age of 10 years at breast-height which corresponds to a plantation age of less than 20 years. This is approximately the termination point of the juvenile corewood or crown-based wood formation process as the plantations achieved crown closure status. These trends and inferences are in general accordance with the observed temporal correlation patterns (e.g., juvenile-mature wood transition occurring at the point of crown closure) and consistent with the developmental patterns reported for other softwood species [24–26].Simultaneously controlling for age-dependent (ontogenetic) and growth-dependent (environmental) influences on fibre attributes is an analytical challenge when quantifying their temporal developmental trends. Attempts to include both ring area and cambial age as predictor variables when specifying fibre attribute prediction models have resulted in mixed success. As shown in Figure1, annual ring area (growth rate) was significantly (p ≤ 0.05) correlated to all 8 attributes for the majority of the years assessed. Consequently, accounting for fluctuations in annual growing conditions on fibre attributes by analyzing temporal developmental patterns by including the chronological calendar year of formation may help in minimizing growth-dependent effects underlying attribute variation. The effect of intraspecific tree competition on growth rate and by extension on fibre attribute variability is partially reduced by stratifying trees into hierarchical size-based classes (e.g., basal area quintiles). Given these inferences, composite multivariate allometric (e.g., [6, 27]) or hierarchical mixed models [28] are plausible analytical frameworks for consideration when attempting to quantify attribute developmental patterns.


### 3.4. Patterns of Attribute Variation and Modeling Implications

As forest management objectives migrate from a volumetric maximization to a value-added product-based paradigm, initiatives to develop fibre attribute prediction models have increased. In regard to coniferous boreal species, the fibre prediction equations developed for jack pine (Pinus banksiana Lamb.) [29], Scot pine (Pinus sylvestris L.) [30]; and black spruce [31] are representative examples of this renewed research focus. Analytically, arriving at the correct specification is largely an iterative process in which observed patterns of attribute development and associated variation across potential covariates are used to select a set of candidate fibre attribute prediction models. In this study, for planted black spruce trees, the variation in fibre attributes was found to decline with increased tree size (Figure 3). Specifically, the coefficient of variation for all 8 attributes decreased from a mean value of approximately 11-12% within the 1st quintile to a mean of 8-9% within the 5th quintile. The range of variation also declined with increasing basal area quintile (tree size). These patterns are suggestive of a size-dependent ordering of attribute variation throughout the stand structure. Combined with the observed size-dependent temporal development trends of the cumulative moving average values, plausible candidate model specifications are those which include size-based explanatory variables. Furthermore, as shown in this study (Figure 1), all 8 attributes were significantly (p≤0.05) correlated with annual ring area (basal area growth rate). Consequently, accounting for fluctuations in annual growing conditions and their corresponding effect on fibre attribute formulation through the inclusion of a growth rate predictor variable may also increase the explanatory power of potential fibre attribute prediction models (e.g., [28]).Figure 3
Fibre attribute variation as measured by the coefficient of variation by size quintile. Note: (1) the solid square denotes the mean value whereas the lower and upper whiskers denote minimum and maximum values, respectively, and wd, ma, me, co, wt, dr, dt, and sa are abbreviative notation for wood density, microfibril angle, modulus of elasticity, fibre coarseness, tracheid wall thickness, radial tracheid diameter, tangential tracheid diameter, and specific surface area, respectively.In summary, the results of this study suggest the fibre attribute models for plantation black spruce should consider the inclusion of predictor variables that reflect age-dependence (intrinsic ontogenetic development (e.g., cambial age)), growth-dependence (climatic-based environmental effects (e.g., annual-ring-width area)), and size-dependence (intraspecific competition (e.g., diameter)) sources of variation.


### 3.5. Relative Quality of Potential Commercial Thinning (CT) Residuals

Given the lack of definitive design specifications for the attributes considered in this study in terms of end-product potential, it was necessary to implement a supplementary comparative analysis in order to attain an appreciation of the quality of the fibre attributes of the black spruce plantations at midrotation. The conceptual rationale for this comparison is that the quality of attributes from more mature and natural-origin black spruce stands would be indicative of the black spruce fibre resource that has historically met or exceeded end-product quality expectations (e.g., [32]). Thus, if attributes from the plantations at midrotation were not substantially different from those observed within natural-origin stands at rotation, it would be reasonable to infer that the end-product potential of residuals arising from plausible commercial thinning treatments would be not substantially dissimilar.Consequently, results from 2 regionally diverse and independent Silviscan-based studies on black spruce attributes were used in this comparison. The first external set was comprised of Silviscan results from 111 breast-height cores obtained from black spruce trees within 75 plots established in 45 stands that were situated on principally upland site-types within Forest Section B4 (Northern Clay) of the Canadian Boreal Forest Region [1] (Source: Dr. B. Pokharel, Nipissing University, North Bay, Ontario, Canada). Although the sample plots were established in stands ranging in age from 25 to 163 years, the majority were within developmentally mature and older stands, as evident from their structural characteristics at time of sampling: mean dominant height of 17.8 m (standard error/minimum/maximum: 3.6/9.0/25.9); mean basal area of 27.4 m2/ha (10.3/17.0/51.4); mean quadratic diameter (cm) of 16.7 (3.8/11.0/28.3); and mean density of 1324 stems/ha (577/1750/2601). The second external data set consisted of results derived from the Silvican analysis of 160 breast-height cores sampled from merchantable-sized (breast-height diameter ≥ 9 cm) black spruce trees located within 16 black spruce stands situated throughout Forest Sections B28a (Grand Falls), B28b (Corner Brook), and B29 (Northern Peninsula) of the Canadian Boreal Forest Region [1] (Source: Table  1 in [33]). Representative mensurational characteristics of the merchantable-sized tree population within these sampled stands suggested that most were at a mature stage of development: mean dominant height of 11.9 m (standard error/minimum/maximum: 2.8/5.9/17.1); mean basal area of 25.2 m2/ha (12.6/1.9/46.3); mean quadratic diameter (cm) of 13.4 (2.5/9.2/18.3); and mean density of 1725 stems/ha (712/250/3200).Denoting these Silviscan data sets from Ontario and Newfoundland as SS-ON and SS-NF, respectively, the results from this comparison revealed slight to moderate differences between the plantations and the mature natural-origin stands in terms of mean lower wood density (487 versus 537 (SS-ON) and 562 (SS-NF)), higher microfibril angle (15.0 versus 12.8 (SS-ON) and 13.9 (SS-NF)), lower modulus of elasticity (13.0 versus 14.9 (SS-ON) and 15.2 (SS-NF)), lower fibre coarseness (332 versus 375 (SS-ON) and 395 (SS-NF)), lower wall thickness (2.3 versus 2.6 (SS-ON) and 2.8 (SS-NF)), higher radial diameter (27.6 versus 27.5 (SS-ON) and 25.9 (SS-NF)), lower tangential diameter (25.9 versus 26.4 (SS-ON) and 26.5 (SS-NF)), and higher specific surface area (336 versus 301 (SS-ON) and 285 (SS-NF)). Similarly, differences among the attributes in terms of their relative variation were not substantially dissimilar: coefficient of variation (%) values of 8.7 versus 9.3 (SS-ON) and 6.8 (SS-NF) for wood density, 24.0 versus 24.3 (SS-ON) and 15.8 (SS-NF) for microfibril angle, 20.8 versus 16.8 (SS-ON) and 11.2 (SS-NF) for modulus of elasticity, 9.4 versus 10.3 (SS-ON) and 7.1 (SS-NF) for fibre coarseness, 8.7 versus 11.5 (SS-ON) and 7.6 (SS-NF) for wall thickness, 7.2 versus 5.1 (SS-ON) and 2.2 (SS-NF) for radial diameter, 5.4 versus 5.3 (SS-ON) and 2.3 (SS-NF) for tangential diameter, and 7.4 versus 9.2 (SS-ON) and 6.7 (SS-NF) for specific surface area.Thus, based on these descriptive comparisons, plantation attributes at midrotation were not substantially different from their mature natural-origin stand counterparts with the possible exception of microfibril angle and fibre coarseness. However, given the importance of these two attributes in terms of their effect on pulp quality (sensu Table1), these differences may be consequential in cases were the thinning residuals are directed towards pulp mills for the production of paper-related end-products. More generally, caution must also be exercised when interpreting such comparative results, in that the size-dependent effect on the developmental patterns of the attributes which was clearly evident for the trees analyzed in this study, is not reflected when employing simple measures of central tendency (sample means).


### 3.6. Potential Management Implications

Although the development of fibre attributes is under strong ontogenetic control and hence largely predetermined, phenotypic plasticity arising from mechanical (support) or hydraulic (physiology) requirements enables trees to modify their developmental patterns while under external stress. This intrinsic plasticity capability helps trees in maximizing their survival probabilities (sensu [34]). As evident by the size-dependent attribute developmental patterns and variability trends observed in this study, black spruce exhibits a considerable range of phenotypic plasticity. Thus, the potential to indirectly influence fibre attribute developmental variation through the regulation of size structures via density control treatments, offers forest managers the opportunity to optimally manipulate not only the quantity but also the quality of the fibre yields. Furthermore, natural pruning, branch diameter growth, juvenile wood production, and the rotational end-products are directly influenced by site occupancy regulation. Specifically, the establishment and maintenance of adequate site occupancies that facilitate rapid crown closure and encourage intraspecific competition have been shown to improve end-product quality at rotation (e.g., black spruce [25], Norway spruce (Picea abies (L.) Karst) [35], Jack pine [36], and Douglas-fir (Pseudotsuga menziesii (Mirb.) Franco) [37]).The size-dependent ordering of fibre attribute development observed by midrotation in this study partially arose from intense intraspecific competition among neighbouring trees once the plantations achieved crown closure status. Asymmetrical competition, principally for above-ground resources (light), has been found to underlie the hierarchical patterns of size-dependent growth within black spruce stands [38]. This competitive effect eventually governs structural development patterns and self-thinning processes. The results of this study suggest that fibre attributes were also adhering to a similar pattern of differentiation as a consequence of competition. Hence, the ability to manipulate competition processes offers forest managers a pathway to indirectly influence fibre attribute developmental patterns and ultimately end-product potentials.


## 3.1. Static Correlative Relationships

The correlation patterns between the fibre attributes largely followed expectation (Table4). Increasing fibre lengths were associated with enhanced lumber related attributes such as lower ma and increased me (stiffness). Decreases in tracheid diameters were associated with increasing wood density and stiffness as evident by the inverse dr–wd-1, dr–me-1, dt–wd-1, and dt–me-1 correlative relationships. Increasing wall thickness which is associated with increased wood density and wood stiffness was also reflected in the results: that is, the positive wt–wd and me–wt correlations. Strong correlations (r≥0.70) were observed among (1) wood density and the modulus of elasticity which are both related to solid wood products and lumber grades (wd–me) and (2) fibre dimensional attributes associated with pulp yields and paper end-products (fl–co, co–wt, co–sa-1, wt–sa-1, and dr–dt). Essentially, these results suggest that increasing wood density is accompanied by a suite of other favourable attributes which could translate into greater solid wood product yields and grades (Table 1). Similarly, long fibre lengths are accompanied with a suite of attributes which collectively are associated with enhanced pulp and paper yields and quality (Table 1).Table 4
Bivariate linear association among fibre attribute values at midrotation as measured by the Pearson product moment correlation coefficient.
Attributea Attributea f l w d m a m e c o w t d r d t f l 1.0 w d 0.0400 1.0 m a −0.3343∗ −0.4840∗ 1.0 m e 0.2995∗ 0.7765∗ −0.8981∗ 1.0 c o 0.6992∗ 0.1969∗ −0.1569 0.2044 1.0 w t 0.4476∗ 0.8201∗ −0.4429∗ 0.6765∗ 0.7201∗ 1.0 d r 0.4799∗ −0.5952∗ 0.2820 −0.4477∗ 0.5913∗ −0.0769 1.0 d t 0.5542∗ −0.4677∗ 0.1585 −0.3037∗ 0.6700∗ 0.0609 0.7072∗ 1.0 s a −0.5090∗ −0.7164∗ 0.3791∗ −0.5812∗ −0.8103∗ −0.9743∗ 0.0560 −0.2047 Note:∗significant correlation at the 0.05 probability level. aDenotations are defined in Table 3.Linear bivariate graphical patterns and significant (p≤0.05) correlations existed between some of the fibre attributes and tree level variables (Table 5). Fibre length was positively correlated with individual-tree diameter, height, and age. Tracheid length is associated with the degree of fibre-bonding and hence is directly related to the tensile and tear strength of derived paper products [20]. Consequently, the observed correlations for fibre length, suggest that the suitability of these plantation black spruce trees for pulp and paper production is enhanced with increasing tree size and age. Similar patterns of correlation between fibre attributes and diameter were observed for height (fl, wd-1, ma, me-1, dr, and dt-H). Positive correlations between age and co, dr, and dt were also detected. Attributes ma and sa were negatively correlated with live crown ratio whereas me was positively correlated. Overall, the attributes were most strongly and frequently correlated the individual-tree size metrics (diameter and height).Table 5
Bivariate linear association between tree level variables and fibre attribute values at midrotation as measured by the Pearson product moment correlation coefficient.
Attributea Tree level variablesb D H A C r f l 0.3636∗ 0.4017∗ 0.4117∗ 0.1808 w d −0.6089∗ −0.4977∗ −0.2559 0.2540 m a 0.3757∗ 0.3337∗ 0.2293 −0.4099∗ m e −0.5144∗ −0.4068∗ −0.2499 0.3752∗ c o 0.2493 0.1898 0.2946∗ 0.1878 w t −0.2832 −0.2330 0.0025 0.2820 d r 0.7186∗ 0.5728∗ 0.4377∗ −0.0671 d t 0.4764∗ 0.3795∗ 0.3477∗ −0.0096 s a 0.1785 0.1704 −0.0296 −0.3185∗ Note:∗significant correlation at the 0.05 probability level. aDenotations are defined in Table 3. bDenotations are defined in Table 2.


## 3.2. Dynamic Correlative Relationships

Figure1 graphically illustrates the pairwise correlation among the attributes on an annual basis (cambial age) as the plantations developed. Six general trends emerged: (1) noncontinuous inverse correlative patterns for the ar–wd, ar–me, ar–wt, wd–dr, me–dr,me–dt,me–sa,and dt–sa relationships; (2) noncontinuous directly proportional patterns for the ar–ma, ar–sa, wd–co, ma–dr, me–wt, co–dr, and co–dt relationships; (3) approximately constant and continuous positive correlations across all cambial ages for the wd–wt and co–wt relationships; (4) approximately constant and continuous negative correlations across all cambial ages for the wd–sa, ma–me, co–sa, and wt–sa relationships; (5) shifting negative-to-positive correlations for the ar–co, ar–dr, ar–dt, wd–me, ma–sa,and dr–dt relationships; and (6) shifting positive-to-negative correlations for the wd–ma, wd–dt, and ma–wt relationships.Figure 1
Matrix of temporal annual correlation patterns for bivariate relationships between individual attributes by cambial age. Note: (1) only Pearson moment correlation coefficients which were significant (p≤0.05) are displayed, and (2) ar, wd, ma, me, co, wt, dr, dt, and sa denote annual-ring-width area, wood density, microfibril angle, modulus of elasticity, fibre coarseness, tracheid wall thickness, radial tracheid diameter, tangential tracheid diameter, and specific surface area, respectively.These results suggest the presence of a set of developmental (age independent) invariant bivariate associations which included thewd–wt, wd–sa, ma–me, co–wt, co–sa, and wt–sa relationships. The underlying morphological (functional) relationships that exist among the attributes may be partially responsible for these results (e.g., wd=fwt). Furthermore, these temporal invariant relationships provide an opportunity to employ surrogate variables for estimating the values of other attributes, through a sequential set of nested bivariate relationships: wt=fwd→sa=fwt→co=fsa.The shifting positive-to-negative and negative-to-positive correlative patterns are plausible transitions points as the plantations reached the crown closure stage of development. This is the stage of stand development that is commonly associated with the cessation of juvenile or crown-based wood formation. These transition points coincided with cambial ages ranging from approximately 8 to 13 years as inferred graphically from Figure1. Assuming a 5-year developmental period is required for the sample trees to reach breast-height (1.3 m), this cambial age range corresponds to a total stand age of approximately 13–18 years. This total age range is similar to the juvenile-mature wood transition age that has previously been reported in the literature for black spruce (i.e., 11–21 year range [21, 22]) and is also similar to that reported for white spruce (Picea glauca (Moench) Voss) [23]. Based on predictions derived from a structural stand density management diagram for these sampled plantations [6], the observed transition age range corresponds to the approximate time of crown closure.Furthermore, the temporal positive-to-negative and negative-to-positive shifts in the correlations may indicate intrinsic differences in the underlying wood formation processes during the juvenile and postjuvenile growth stages. The majority of the attributes which exhibited transitions in their correlative patterns were directly or indirectly associated with cell geometry; for example, increasing radial diameters were associated with decreasing tangential diameters during the juvenile phase whereas they were directly proportional to each other during the maturing phase. Results derived from other studies are in accordance with this inference. For example, Lenz et al. [23] suggested that the accelerated growth which is likely the outcome of frequent anticline divisions of cambial initials partially underlies the unique geometrical characteristics of cells produced during the juvenile phase.Collectively, these results provide confirmatory support for previous research in which cambial age has been found to be associated with unique fibre formation patterns [24]. For example, microfibril angle has been used as an investigative tool to differentiate between juvenile and mature wood types [25]. The higher growth rate and production of juvenile wood during the precrown closure stage may also elicit a temporary decline in wood quality in terms of end-product potential (e.g., lower wd and me and higher ma (Table 1)). Overall, the results and supporting references suggest that physiologically induced changes in allometric relationships were partially responsible for the observed shift in correlation patterns as the black spruce plantations transitioned from the precrown closure stage of development to the postcrown closure stage of development.


## 3.3. Temporal Developmental Patterns

The temporal development trends of the cumulative area-weighted moving average values by tree size (basal area quintile) for wood density, microfibril angle, modulus of elasticity, coarseness, tracheid wall thickness, tracheid radial diameter, tracheid tangential diameter, and specific surface area are graphically presented in Figure2. Irrespective of attribute, the patterns revealed that the most rapid period of change was observed within the first decade of growth and that the rate of change varied inversely with tree size (e.g., rate of change declined with increased tree size). Among the attributes, Wd, Co, Wt, Dr, Dt, and Sa exhibited the greatest rate of change during this period whereas Ma and Me exhibited the least. Following this initial period of rapid change, (1) Wd, Co, Wt, Dr, Dt, and Sa exhibited a period of relative stability that was characterized by a marginal but continuous increase (Wd, Co, Wt, Dr, and Dt) or decrease (Sa) with increasing cambial age, and (2) Ma and Me exhibited a continuous declining (Ma) or increasing (Me) trend with increasing cambial age. For a given cambial age, Wd, Me, Co, and Wt exhibited size-dependency in which values declined with increasing basal area quintile. After approximately 40 years, the rate of change within a given basal area quintile was minimal.Figure 2
Three-dimensional visualization of mean temporal developmental trends illustrating patterns of change in attribute-specific cumulative area-weighted moving average values with increasing cambial age and size quintile. Note:Wd, Ma, Me, Co, Wt, Dr, Dt, and Sa denote the mean value up to the ith cambial age for wood density, microfibril angle, modulus of elasticity, fibre coarseness, tracheid wall thickness, radial tracheid diameter, tangential tracheid diameter, and specific surface area, respectively.The results indicated that the most rapid change within the trajectories occurred before the trees reached a cambial age of 10 years at breast-height which corresponds to a plantation age of less than 20 years. This is approximately the termination point of the juvenile corewood or crown-based wood formation process as the plantations achieved crown closure status. These trends and inferences are in general accordance with the observed temporal correlation patterns (e.g., juvenile-mature wood transition occurring at the point of crown closure) and consistent with the developmental patterns reported for other softwood species [24–26].Simultaneously controlling for age-dependent (ontogenetic) and growth-dependent (environmental) influences on fibre attributes is an analytical challenge when quantifying their temporal developmental trends. Attempts to include both ring area and cambial age as predictor variables when specifying fibre attribute prediction models have resulted in mixed success. As shown in Figure1, annual ring area (growth rate) was significantly (p ≤ 0.05) correlated to all 8 attributes for the majority of the years assessed. Consequently, accounting for fluctuations in annual growing conditions on fibre attributes by analyzing temporal developmental patterns by including the chronological calendar year of formation may help in minimizing growth-dependent effects underlying attribute variation. The effect of intraspecific tree competition on growth rate and by extension on fibre attribute variability is partially reduced by stratifying trees into hierarchical size-based classes (e.g., basal area quintiles). Given these inferences, composite multivariate allometric (e.g., [6, 27]) or hierarchical mixed models [28] are plausible analytical frameworks for consideration when attempting to quantify attribute developmental patterns.


## 3.4. Patterns of Attribute Variation and Modeling Implications

As forest management objectives migrate from a volumetric maximization to a value-added product-based paradigm, initiatives to develop fibre attribute prediction models have increased. In regard to coniferous boreal species, the fibre prediction equations developed for jack pine (Pinus banksiana Lamb.) [29], Scot pine (Pinus sylvestris L.) [30]; and black spruce [31] are representative examples of this renewed research focus. Analytically, arriving at the correct specification is largely an iterative process in which observed patterns of attribute development and associated variation across potential covariates are used to select a set of candidate fibre attribute prediction models. In this study, for planted black spruce trees, the variation in fibre attributes was found to decline with increased tree size (Figure 3). Specifically, the coefficient of variation for all 8 attributes decreased from a mean value of approximately 11-12% within the 1st quintile to a mean of 8-9% within the 5th quintile. The range of variation also declined with increasing basal area quintile (tree size). These patterns are suggestive of a size-dependent ordering of attribute variation throughout the stand structure. Combined with the observed size-dependent temporal development trends of the cumulative moving average values, plausible candidate model specifications are those which include size-based explanatory variables. Furthermore, as shown in this study (Figure 1), all 8 attributes were significantly (p≤0.05) correlated with annual ring area (basal area growth rate). Consequently, accounting for fluctuations in annual growing conditions and their corresponding effect on fibre attribute formulation through the inclusion of a growth rate predictor variable may also increase the explanatory power of potential fibre attribute prediction models (e.g., [28]).Figure 3
Fibre attribute variation as measured by the coefficient of variation by size quintile. Note: (1) the solid square denotes the mean value whereas the lower and upper whiskers denote minimum and maximum values, respectively, and wd, ma, me, co, wt, dr, dt, and sa are abbreviative notation for wood density, microfibril angle, modulus of elasticity, fibre coarseness, tracheid wall thickness, radial tracheid diameter, tangential tracheid diameter, and specific surface area, respectively.In summary, the results of this study suggest the fibre attribute models for plantation black spruce should consider the inclusion of predictor variables that reflect age-dependence (intrinsic ontogenetic development (e.g., cambial age)), growth-dependence (climatic-based environmental effects (e.g., annual-ring-width area)), and size-dependence (intraspecific competition (e.g., diameter)) sources of variation.


## 3.5. Relative Quality of Potential Commercial Thinning (CT) Residuals

Given the lack of definitive design specifications for the attributes considered in this study in terms of end-product potential, it was necessary to implement a supplementary comparative analysis in order to attain an appreciation of the quality of the fibre attributes of the black spruce plantations at midrotation. The conceptual rationale for this comparison is that the quality of attributes from more mature and natural-origin black spruce stands would be indicative of the black spruce fibre resource that has historically met or exceeded end-product quality expectations (e.g., [32]). Thus, if attributes from the plantations at midrotation were not substantially different from those observed within natural-origin stands at rotation, it would be reasonable to infer that the end-product potential of residuals arising from plausible commercial thinning treatments would be not substantially dissimilar.Consequently, results from 2 regionally diverse and independent Silviscan-based studies on black spruce attributes were used in this comparison. The first external set was comprised of Silviscan results from 111 breast-height cores obtained from black spruce trees within 75 plots established in 45 stands that were situated on principally upland site-types within Forest Section B4 (Northern Clay) of the Canadian Boreal Forest Region [1] (Source: Dr. B. Pokharel, Nipissing University, North Bay, Ontario, Canada). Although the sample plots were established in stands ranging in age from 25 to 163 years, the majority were within developmentally mature and older stands, as evident from their structural characteristics at time of sampling: mean dominant height of 17.8 m (standard error/minimum/maximum: 3.6/9.0/25.9); mean basal area of 27.4 m2/ha (10.3/17.0/51.4); mean quadratic diameter (cm) of 16.7 (3.8/11.0/28.3); and mean density of 1324 stems/ha (577/1750/2601). The second external data set consisted of results derived from the Silvican analysis of 160 breast-height cores sampled from merchantable-sized (breast-height diameter ≥ 9 cm) black spruce trees located within 16 black spruce stands situated throughout Forest Sections B28a (Grand Falls), B28b (Corner Brook), and B29 (Northern Peninsula) of the Canadian Boreal Forest Region [1] (Source: Table  1 in [33]). Representative mensurational characteristics of the merchantable-sized tree population within these sampled stands suggested that most were at a mature stage of development: mean dominant height of 11.9 m (standard error/minimum/maximum: 2.8/5.9/17.1); mean basal area of 25.2 m2/ha (12.6/1.9/46.3); mean quadratic diameter (cm) of 13.4 (2.5/9.2/18.3); and mean density of 1725 stems/ha (712/250/3200).Denoting these Silviscan data sets from Ontario and Newfoundland as SS-ON and SS-NF, respectively, the results from this comparison revealed slight to moderate differences between the plantations and the mature natural-origin stands in terms of mean lower wood density (487 versus 537 (SS-ON) and 562 (SS-NF)), higher microfibril angle (15.0 versus 12.8 (SS-ON) and 13.9 (SS-NF)), lower modulus of elasticity (13.0 versus 14.9 (SS-ON) and 15.2 (SS-NF)), lower fibre coarseness (332 versus 375 (SS-ON) and 395 (SS-NF)), lower wall thickness (2.3 versus 2.6 (SS-ON) and 2.8 (SS-NF)), higher radial diameter (27.6 versus 27.5 (SS-ON) and 25.9 (SS-NF)), lower tangential diameter (25.9 versus 26.4 (SS-ON) and 26.5 (SS-NF)), and higher specific surface area (336 versus 301 (SS-ON) and 285 (SS-NF)). Similarly, differences among the attributes in terms of their relative variation were not substantially dissimilar: coefficient of variation (%) values of 8.7 versus 9.3 (SS-ON) and 6.8 (SS-NF) for wood density, 24.0 versus 24.3 (SS-ON) and 15.8 (SS-NF) for microfibril angle, 20.8 versus 16.8 (SS-ON) and 11.2 (SS-NF) for modulus of elasticity, 9.4 versus 10.3 (SS-ON) and 7.1 (SS-NF) for fibre coarseness, 8.7 versus 11.5 (SS-ON) and 7.6 (SS-NF) for wall thickness, 7.2 versus 5.1 (SS-ON) and 2.2 (SS-NF) for radial diameter, 5.4 versus 5.3 (SS-ON) and 2.3 (SS-NF) for tangential diameter, and 7.4 versus 9.2 (SS-ON) and 6.7 (SS-NF) for specific surface area.Thus, based on these descriptive comparisons, plantation attributes at midrotation were not substantially different from their mature natural-origin stand counterparts with the possible exception of microfibril angle and fibre coarseness. However, given the importance of these two attributes in terms of their effect on pulp quality (sensu Table1), these differences may be consequential in cases were the thinning residuals are directed towards pulp mills for the production of paper-related end-products. More generally, caution must also be exercised when interpreting such comparative results, in that the size-dependent effect on the developmental patterns of the attributes which was clearly evident for the trees analyzed in this study, is not reflected when employing simple measures of central tendency (sample means).


## 3.6. Potential Management Implications

Although the development of fibre attributes is under strong ontogenetic control and hence largely predetermined, phenotypic plasticity arising from mechanical (support) or hydraulic (physiology) requirements enables trees to modify their developmental patterns while under external stress. This intrinsic plasticity capability helps trees in maximizing their survival probabilities (sensu [34]). As evident by the size-dependent attribute developmental patterns and variability trends observed in this study, black spruce exhibits a considerable range of phenotypic plasticity. Thus, the potential to indirectly influence fibre attribute developmental variation through the regulation of size structures via density control treatments, offers forest managers the opportunity to optimally manipulate not only the quantity but also the quality of the fibre yields. Furthermore, natural pruning, branch diameter growth, juvenile wood production, and the rotational end-products are directly influenced by site occupancy regulation. Specifically, the establishment and maintenance of adequate site occupancies that facilitate rapid crown closure and encourage intraspecific competition have been shown to improve end-product quality at rotation (e.g., black spruce [25], Norway spruce (Picea abies (L.) Karst) [35], Jack pine [36], and Douglas-fir (Pseudotsuga menziesii (Mirb.) Franco) [37]).The size-dependent ordering of fibre attribute development observed by midrotation in this study partially arose from intense intraspecific competition among neighbouring trees once the plantations achieved crown closure status. Asymmetrical competition, principally for above-ground resources (light), has been found to underlie the hierarchical patterns of size-dependent growth within black spruce stands [38]. This competitive effect eventually governs structural development patterns and self-thinning processes. The results of this study suggest that fibre attributes were also adhering to a similar pattern of differentiation as a consequence of competition. Hence, the ability to manipulate competition processes offers forest managers a pathway to indirectly influence fibre attribute developmental patterns and ultimately end-product potentials.


## 4. Conclusions

The objectives of this study were to analyze the temporal developmental patterns of commercially relevant fibre attributes and their interrelationships within maturing black spruce plantations. The results from the analyses indicated that fibre attributes were correlated among themselves and with morphological tree characteristics and exhibited unique size-dependent temporal developmental patterns. After approximately 40 years of development, the rate of change for the majority of the attributes was minimal. The period of most rapid change was during the precrown closure period of stand development which encompassed the first 20 years of growth. A pattern of size-dependent ordering of attribute variation within the black spruce plantations was evident from the analyses. Although the results of this study suggested that commercial row thinning of similar density-stressed plantations at midrotation may potentially yield satisfactory fibre quality for specific end-products,the lack of product-based design specifications negated a conclusive determination. Identification of some of the sources of variation underlying the development of the studied attributes enabled inferences to be derived for potential use in future modeling studies.



---

*Source: 7895289-2016-11-08.xml*