Assessing Soil-Related Black Spruce and White Spruce Plantation Productivity

This article focuses on modelling and mapping the productivity of black (Picea mariana) and white spruce (Picea glauca) plantations across the Black Brook forest management area in northwestern New Brunswick, Canada, encompassing about 200,000 ha. This effort involved establishing 3500 50 m survey plots, each informing about: plantation age (15 to 43 years), planted species type, stem count, tree height, basal area, and wood volume. All of this was supplemented with location-specific productivity predictors, i.e., xy location and specifications pertaining to soil type, soil drainage (established through digital elevation modelling by way of the depth-to-water index DTW), and years since thinning (pre-commercial and commercial), and. The DTW index, as it emulates the elevation rise away from open water features such as streams, rivers and lakes, allowed the re-mapping of existing soil borders by topographic position and drainage association. Non-linear regression analysis revealed that plantation height, basal area and volume all increased with plantation age, as to be expected. Pre-commercial thinning in plantations <30 years old had a positive while the more recent commercial thinning still had the negative effect on standing wood volume and mean annual volume increment (MAI). White spruce MAI generally exceeded black spruce (MAI) by a factor of 1.25. Poor and excessive soil drainage reduced MAI. Best growth performances occurred on plantations established on well-drained calcareous soils. The best-fitted results so obtained allowed for generating black and white spruce MAI maps across the forest management area by ridge-to-valley soil and DTW location at 10 m resolution. These maps were subsequently used for site-by-site silvicultural evaluation and ranking purposes.


Introduction
In general, of forest plantation productivities are not only affected by forest management actions including species selection and stocking level to promote plantation growth and yields but also by climate and soil conditions.This article focuses on aspects pertaining to how differences in soil quality and drainage affect the growth of black and white spruce plantations, with soil type referring to differences in landforms and related surface expressions and geological substrates.Likely soil-related growth and yield affecting differences refer to rooting depth, organic matter and coarse fragment content, and nutrient availabilities as modified by rate of soil weathering, litter decomposition, soil moisture regime and climate.Climate related factors refer to length of growing season, the extent of atmospheric deposition in the form of precipitation (rain, snow, and nutrients such as N, S, Ca, Mg and K), and actual evapotranspiration as dictated by monthly precipitation and air temperature patterns (Keys et al., 2016).
Highest soil productivities can be expected to occur on deep, friable and organically enriched soils with well-functioning root systems.In other soils, root growth is often restricted by high coarse fragment content, soil compaction, lack of aeration, low temperature, and/or the presence of toxic substances (Table 1, Arp & Krause, 2002).Typically, plantation growth is better maintained by organically enriched soils on well-drained locations due to prolonged retention of plant available water.On poorly to imperfect drainage locations, insufficient soil aeration due to persisting water saturation restricts root development (Kozlowski, 1985;Gale and Grigal, 1987;Drever & Lertzman, 2001;Crow, 2005;Angstmann et al., 2013).For example, Douglas fir and white spruce do not tolerate low soil aeration, but black spruce, red alder, western cedar, Sitka spruce and western hemlock are adapted to shallow water tables.Rooting is also affected Table 1.General relationships between soil properties and site quality, including optimal conditions.by soil temperatures (Alvarez-Uria and Körner, 2007), being slow and delayed in soils that remain cool as seasons transit from winter to spring and summer.In addition, fine roots are sensitive to acidification, mainly due to acid-induced Al mobilization (Cronan and Grigal, 1994).
This article focuses on stem count, tree height, basal area and wood volume growth in select black and white spruce plantations across a 250,000 ha forest management area in northwestern New Brunswick (Figure 1).The research involved discerning how these growth metrics are affected by species type (white versus black spruce), soil type and soil drainage.Soil type was defined by type-specific landform and geological substrate combinations, also referred to as soil associations or catenas.Soil drainage type was adjusted to conform to the DTM-derived depth-to-water pattern (DTW) across the land as generated from a 10-m resampled LiDAR-generated digital terrain model (DTM) for the area.

Study Area
The forest management area in northwestern New Brunswick refers to the Black Brook management zone (longitude −67.66; latitude 47.231) spans elevations from 128 to 569 above sea level (Figure 1).Mean annual precipitation amounts to 1104 mm; mean annual air temperature is 3.5˚C; mean degree growing days above 5˚C amount to 1532.6 (Environment Canada, 2017).The dominant upland species are Sugar Maple (Acer saccharum) (30%), Yellow Birch (Betula alleghaniensis) (11%), and Balsam Fir (Abies balsamea) (10%).Planted stands cover 33% of the management zone and consist off White Spruce (Picea glauca, 42%), Black Spruce (Picea mariana, 40%), and Norway Spruce (Picea abies, 14%).Surface deposits vary from alluvium within floodplains, glacio-fluvial along the floodplains, residuals on ridges and steep slopes, and combinations of ablation and/or basal till along generally undulating to rolling terrain (Langmaid et al. 1976;Colpitts et al., 1995;Fahmy et al., 2002).An overview of forest soil units and associated bedrock formations is provided in Table 2.

DTM-Based DTW Delineations
The DTM raster at 1-m resolution was generated from the 2014 LiDAR pointcloud data for the area, through systematic point-of-last-return classification, followed by natural neighbour interpolation.This DTM was then re-sampled at 10 m resolution and used to generate the cartographic depth-to-water delineation as developed from the DTM-generated datalayers (Murphy et al., 2011).

Soil Association and Mode of Surface Deposition: Shapefile Adjustments
The DTM-derived slopes and DTW patterns away from the stream and floodplain flow channels were used to hydro-topographically correct the existing soil catena borders to conform to the lake, river, and wetland features of the area.In

Plantation Survey
In 2007, 3500 50-m 2 circular forest development survey plots were used to generate data for stem count, tree height, basal area and stem volume for select black and white spruce plantations across the Black Brook forest management area.These plots were systematically located within select plantation along a 100 × 100 m sampling grid.All planted and non-planted trees were tallied for each plot by species, diameter at breast height.Plot height was set to be the height of tree closest to each grid centre.These determinations were used to estimate stem density (stems/ha), basal area (m 2 /ha), total volume (m 3 /ha) and mean annual increment (m 3 /ha/year) for each plot.Also noted was plantation age and years since pre-commercial (n = 791) and commercial (n = 203) thinning.Each plot was further characterized by planted species, landform by mode of deposition, soil type, elevation, slope, aspect, and DTW.The number of plots per soil type are listed in Table 2. Due to the nature of the plot survey and the area limitations per soil type, observations per soil type varied widely, from 6 to 1291.In addition, white spruce plantations (n = 1050) were generally established on well drained sites.

Data Management
All plot-based numerical and categorical data were compiled into a single spreadsheet, with column headings identifying all the above variables, and each row providing the data entries, one plot at a time.Planted species and soil type per plot were binary coded (present 1; not present 0).The plot-specific numbers for DTW numbers were transformed into log 10 DTW to capture the changing drainage conditions from very wet to dry in a linear fashion.Also added where the DEM-extracted values for slope and elevation for each plot at centre.

Regression Model and Analysis
The compiled data were subject to multivariate nonlinear regression analysis (stepwise backward).Plot-estimated per hectare stem counts, height, basal area, total volume and mean annual total volume increment (MAI = total wood volume / planation age) were used as dependent variables.All the other numerical and binary variables were entered as independent variables.The regression model was formulated as follows: for stems: for tree height, basal area and total wood volume: where WS refers to white spruce (coded 1 when present, zero otherwise for black spruce), ( ) ( ) ( ) a, b, c e, d, f, g, h, k are adjustable coefficients, CL and CT referring to plantations with pre-commercial and commercial thinning operations, and SF as the soil factor, given by with soil types numbered soil1, soil2, soil3 , coded 1 or 0 for presence and absence, respectively, and with 1 2 3 , , s s s  as adjustable parameters.
Equations ((1a) and ( 1b)) differ because stem counts tend to decrease with increasing age due to self-thinning and other processes, while tree height, basal area and total stem volume increase asymptotically with age towards upper species-, climate-and soil-moderated limits.The step-wise regression analysis systematically eliminated variables with non-significant predictor coefficients due to low probability (P < 0.001) and t-values (t < 3).The best-fitted regression models served to generate the predictor equations for stem count, tree height, basal area and total tree volume.

Plantation Productivity Mapping, by Plantation Type
The best-fitted model for MAI was used to map MAI across the entire forest management area, using ArcMap procedures guided as follows: 1) Two maps were generated, one for black spruce and one for white spruce.
2) Plantation age was fixed at 20 years.
3) Stem count was fixed at 1 2500 stems ha − ⋅ , to reflect typical planting stock density.
4) The best-fitted soil coefficients were assigned to each soil type (catena) where present.
5) The effects of pre-commercial and commercial thinning were not considered.

Analytical Limitations
Even with the availability of the 3500 plot determinations, there is the issue of non-equal representation by plantation type (n = 1050 for black spruce; n = 2450 for white spruce) and soil type (Table 2).In addition, there is an uneven plot distribution across the DTW classes, with the majority located within the 1 < DTW < 10 m.As result, it was not possible to determine how plantation growth would be affected by species-specific sensitivities to soil type and topographic position outside this DTW range.For that reason, the h and k parameters in Equation (2) were set equal to h = 2.6, and k = 4.8 for both species, to enable growth predictions across the entire DTW range from essentially open water conditions to steep ridge tops.In addition, all soil parameters were kept common not only by species but also for stem count, tree height, basal area and total wood volume.

Results
In total, there were 2450 black spruce and 1050 white spruce plots (Figure 1).
The plot-generated variables are presented in Table 3 by listing the basic statistics for plantation age, height, basal area, total volume, MAI, and for the DEM-derived variables log 10 DTW, slope, elevation.The histograms of these variables (Figure 2) show wide variations for each these variables, and especially so for the stem-counts: from 20 to 20,400 per hectare.The higher stem count numbers were due to dense non-planted species regenerations.Low stem counts occurred on wet spots, on shallow soils with rock exposures, and densely machine-trafficked areas.A similar stem-count range was reported by Pitt and Lanteigne (2008) in their work on un-thinned and thinned balsam fir and red

Best-Fitted Regression Coefficients
The best-fitted regression coefficients for stem count tree height, basal area and total volume by way of Equation ( 1) are listed in Table 4 and Table 5.Also shown are their least-squares coefficient estimates, and associated errors of esti- tions.This would be due to the intended basal area and wood volume reductions at the time of commercial thinning, and could in part also be due to post-commercial thinning mortality (Thorpe et al., 2008).
With respect to soil drainage, stem count, tree height, basal area and total wood volume were all similarly affected, being optimal within the 1 < DTW < 10 m range, with tree height, basal area and wood volume productivity somewhat more sensitive to extreme soil moisture regimes than stem count.This is in agreement with Wang and Klinka (1996)  With respect to soil type, there were also significant differences for stem count, tree height, basal area and total wood volume, as listed in Table 5.In particular, this list revealed that the effects of soil type on wood volume productivity were most strongly correlated to the soil type contributions pertaining to basal area and tree height, but were only weakly correlated to the soil contributions pertaining to stem count (Table 6).Sorting the soil contributions to wood volume (or MAI) from low to high indicated that the greatest gains in wood volume occurred on well-drained and medium texture calcareous soils (Figure 3).The lowest gains occurred on coarse-textured and gravelly outwash plain soils (Muniac) and fine-textured basal till soils (Siegas).Among the soils on siliceous substrate, plantation productivity decreased with increasing stoniness: least for Holmesville, and most pronounced for Victoria.The ranking of the best-fitted regression coefficients between total wood volume and soil type generally agrees with the soil property descriptions in Table 2, whereby the Caribou, Undine and Kedgwick soils encompass most if not all of the optimal soil conditions as per Table 1.The opposite applies to the Muniac and Siegas soils.
The multivariate regression results varied in the amount of plot-to-plot variance explained, i.e., R 2 = 0.228 for stem counts, 0.588 for tree height, 0.399 for basal area, and 0.530 for tree volume (Table 4).Hence, much of the plot-to-plot variations remain unresolved by way of Equations ( 1)-(3), and especially so for stem count.Choosing MAI as dependent variable explained 33% of the overall MAI variations.The spread of the regression-generated best-fitted residuals for tree height, basal area, total wood volume and MAI are plotted in Figure 4 and Table 6.Best-fitted multiple regression results for the best-fitted soil type contributions to total wood volume (Table 5), using the best-fitted soil type contributions for basal area, tree height and stem count (Table 5) as predictor variables.Figure 3. Arranging the regression coefficients between total wood volume and soil type (Table 5) from low to high.for MAI (Figure 6).

Limitations
The above results do not match the precision that is generally needed for forest related growth and yield predictions.This is in part due to: i) Procedural issues, whereby the plots were selected for a one-time grid evaluation of the overall plantation growth response variations.For best growth and yield results, plots need to be re-measured at pre-determined intervals.
ii) The age-limited distribution of the plantations.This limitation does not al- iii) The unequal distribution of the plots by age, species and soil type and drainage conditions.This limits the extent to which the above results fully represent the entire variations across these conditions.For example, the results obtained    for the less frequent soil types would be more influenced by chance occurrences and outliers than the results obtained for the more frequent soil types.iv) The need to map growth-affecting soil properties, from ridge tops to valleys mapping, with soil properties referring to soil rooting depth, texture, organic matter, coarse fragment content and mineralogy.Doing this requires further in-depth studies and model calibrations.v) For generalization purposes, the evaluations need to be placed into the context of regionally varying climate and atmospheric deposition conditions (Keys et al., 2016).
In addition to the above plot-based evaluations, one could use LiDAR-generated forest metrics for stem counts, tree height, basal and total wood volume at the tree-by-tree level, with, e.g., 20 × 20 or 25 × 25 m 2 summations (Strunk et al. 2008;Woods et al., 2011;Treitz et al., 2012;Wulder et al., 2012;Latifi et al., 2015).In doing this, Tompalski et al., 2016 reported the following relative root mean square differences between air-borne laser scanned and stand-matched yield curve projections regarding select forest stands on Vancouver Island: for maximum tree height, 31.1%; for basal area, 19.8%; for total stem volume, 21.8%; all at a stand-level age of 80 years.This was done through single-tree summations for 9559.25 × 25 m 2 grid plots representing four western species groups (hemlock, cedar, fir/spruce, alder).Overall tree height and stem wood volumes amounted to 33 m, and 576 m 3 •ha −1 .
Since productivity may range from low to high on imperfectly to poorly drained soil based on micro-topographic mound to pit variations, it would be important to examine this aspect as well.However, whether this detail can be generated from LiDAR-DTMs remains to be seen, because point densities for LiDAR-generated ground reflections can be quite low underneath densely spaced forest plantations.It may, however be possible to infer the extent of mound-and pit surface expressions based on soil type characteristics that lead to wind-induced tree uprooting due to shallow rooting.

MAI Productivity Mapping
The best-fitted Equations ( 1)-(3) were used to produce white and black spruce MAI productivity maps at age 20 across the entire Black Brook forest management area, as shown in Figure 7, with excerpts in Figure 8(a) and Figure 8(b).Overlaid on these maps are the plot-based black and white spruce plantation MAI determinations.In general, these maps illustrate that the overall landscape dependence of MAI is-apart from plantation age and species-mainly affected by soil type and drainage.In detail, MAI is lowest in wet spots, is best along the lower portions of rolling to hummocky slopes, and gradually decreases from there towards steeper ridge tops.The positive effect of soil type on MAI is most pronounced on deep and well-drained soils with calcareous content, such as the Caribou, Jardine, Kedgwick and Undine soils.

Concluding Remarks
Examining the results of the 3500 forest development survey plots in relation to Figure 8. Actual versus mapped MAI conformance details from Figure 7 for the black and and white spruce plantations (left and right sides of panels, respectively) within the Black Brook Forest Management area.Top: part of the top eastern management area, centered on mostly calcareous soils (Caribou, Jardine).Bottom: part of the mid-western management area, centered on siliceous soils (Glassville, Holmesville, Victoria, Grand Falls).plantation age, planting stock, thinning operations, soil type and drainage allowed for the development of a generalizable model (Equations ( 1)-( 3)).The resulting best-fitted model was used to estimate and predict the productivities of black and white spruce plantations by soil type, from ridge tops to valleys based on theoretical considerations: e.g., no growth on wet soils and fully exposed bedrock, and varying growth expectations by soil type as affected by substrate lithology, texture and accessible rooting depth.In view of currently forthcoming forest metrics that detail tree-by-tree performances in terms of stem count, height, basal area and wood volume within, e.g., each 20 × 20 m 2 square across entire forest management areas, more of this could be done at a much finer resolution.Doing so is promising because-even in view of the above limitations-about one third to one half of the plot-based productivity variations for tree height, basal area and total wood volume can be quantified using Equations ( 1)-(3).While the results so generated are still insufficient for precise growth and yield predictions across the plantations and beyond the existing age range, they are providing a starting point for further soil factor quantifications in forest management.At this stage, the mapped white and black spruce projections in Figure 6 and Figure 7 provide a semi-quantitative bases for stand-and sitespecific ranking purposes, and-as such-are being used for directing silvicultural operations (harvesting, site preparations, seedling stock selection planting, thinning) away from low productivity sites that are either too wet or too dry.

Figure 1 .
Figure 1.Black Brook forest management area within northwestern New Brunswick, Canada: location, extent, terrain elevation, black and white spruce plantation survey plots, and distribution of soils with calcareous content.
This delineation involved the following raster-based algorithms: fill, flow direction, flow accumulation, stream network classification (4 ha for stream flow initiation), and determining the cartographic depth-to-water (DTW) pattern as dictated by the least elevation rise away from the stream network pattern for which DTW = 0. Floodplains were delineated in the same way, by setting the upslope flow-initiation area for floodplain initiation requiring about 400 ha of upstream flow accumulation areas.Each floodplain so located was subsequently delineated by limiting the floodplain elevation rise to DTW = 4 m away from the DTW = 0 defining floodplain channels.
so doing, open-water features were excluded from the soil polygons, and organic soils would be relocated and/or limited to the extent of the corresponding image-recognized wetland features.Residual soil formations were constrained to occur along steep ridges and slopes only.The lower-lying borders of the till and glacio-fluvial formations were set to coincide with their adjacent DTM-recognized floodplains borders.The existing soil drainage specifications within each soil catena were replaced by the stream-channel based DTW-pattern using the following specifications: <10 cm, very poor; 10 to 25 cm, poor; 25 to 50 cm, imperfect; 50 to 100 cm, moderate; 1 to 20 m well, >20 m excessive.In addition, the geological bedrock formation outlines were adjusted to conform to the lithology descriptions per soil association polygon.The soil types so adjusted are listed in
who provided the metric tree height index assessment of Ontario in terms of the soil moisture regime variations: moderately dry, 15.8 m; slightly dry, 20 m; fresh, 20.5 m; moist, 20 m; very moist 16.3 m; wet, 13 m; very wet, 5.8 m.

Figure 5
Figure 5 by way of scatterplots and histograms.These results are similar to what is normally found across plot-based modelling efforts, as documented by, e.g.,Pinjuv et al. (2007) in Dykstra(2007).Using actual instead of modelled stem counts improved the analytical results slightly by increasing regression-captured variations towards R 2 values amounting to 0.613 for basal area, 0.676 for total wood volume, and 0.434 for MAI.This small increase is not surprising since the plot-based basal area and volume estimates already include the plot-generated stem-count estimates.By plot-to-model comparisons, the following 80% conformance results were obtained: ±1300 stems ha −1 , ±3.8 m for tree height, ±11 m•ha −1 for basal area, ±55 m 3 •ha −1 for total wood volume, and 2.2 m 3 •ha −1 •year −1 therefore set at 3 to allow for the curvilinear fit for height, basal area and total wood volume within the 0 to 40-year range.The extent to which this determination remains valid beyond age 40 requires further research.

Figure 4 .
Figure 4. Best-fitted scatterplots and histograms for the best-fitted residuals for basal area and total wood volume (n = 3500).

Figure 5 .
Figure 5. Best-fitted scatterplots and histograms of the best-fitted residuals for plantation height, and MAI (n = 3500).

Figure 6 .
Figure 6.Conformance plots showing the cumulative frequencies of the absolute residual between the individual plot determinations and corresponding best-fitted model results for stem count, tree height, basal area, total wood volume, and mean annual increment (MAI).

Table 2 .
Overview oFahmy et al., 2010)associated bedrock geology, mode of deposition, and lithology of parent material for the Black Brook management zone, northwestern New Brunswick, Canada.Bedrock Geology retrieved from (New Brunswick Department of Natural Resources N.B.DNR, 2008;Fahmy et al., 2010).

Table 2 ,
by mode of surface deposition, texture, and by lithology with emphasis on calcareous versus siliceous content.

Table 3 .
Variables informing about site-specific attributes based on plot observations (n = 3500).

Table 4 .
Best-fitted regression coefficients for stem count tree height, total wood volume, and basal area by plantation age, white versus black spruce, soil moisture regime (quantified by way of DTW), and thinning intervention (all p values <0.001 except larger where absolute t-values <4).

Table 5 .
Best-fitted regression coefficients for stem count tree height, total wood volume, and basal area by soil type (all p-values < 0.001).