Modelling Climate Effects on Site Productivity and Developing Site Index Conversion Equations for Jack Pine and Trembling Aspen Mixed Stands

Abstract

Forest site productivity estimates are crucial for making informed forest resource management decisions. These estimates are valuable both for the tree species currently growing in the stands and for those being considered for future stands. Current models are generally designed for pure stands and do not account for the influence of climate on tree growth. Consequently, site index (SI) conversion equations were developed specifically for jack pine (Pinus banksiana Lamb.) and trembling aspen (Populus tremuloides Michx.) trees grown in naturally originated mixed stands. This work involved sampling 186 trees (93 of each species) from 31 even-aged mixed stands (3 trees per species per site) across Ontario, Canada. Stem analysis data from these trees were utilized to develop stand height growth models by incorporating climate variables for each species. The models were developed using a mixed effects modelling approach. The SI of one species was correlated with that of the other species and climate variables to establish SI conversion equations. The effect of climate on site productivity was evaluated by projecting stand heights at four geographic locations (east, center, west, and far west) in Ontario from 2022 to 2100 using the derived stand height growth models. Height projections were made under three emissions scenarios reflecting varying levels of radiative forcing by the end of the century (2.6, 4.5, and 8.5 watts m⁻²). Climate effects were observed to vary across different regions, with the least and most pronounced effects noted in the central and far western areas, respectively, for jack pine, while effects were relatively similar across all locations for trembling aspen. Stand heights and SIs of jack pine and trembling aspen trees grown in naturally originated mixed stands can be estimated using the height growth models developed here. Similarly, SI conversion equations enable the estimation of the SI for one species based on the SI of another species and environmental variables.

1. Introduction

Site productivity is the composite expression of a variety of physical and chemical attributes of a forested area, including its soil, topography, and climate. It affects tree growth, recruitment, and mortality in a stand [1]. Therefore, forest site productivity estimates are very important in making informed forest management decisions. These estimates are not only important for the tree species currently growing in the stands but also for the ones considered to be growing in the future. The site index (SI) defined as the stand height (mean height of codominant and dominant trees with their crowns extending at or above the general level of canopy) at a specified stand age is a commonly used measure of site productivity [2]. It is the main driver of most growth and yield models that are used to estimate wood supply and prepare forest management plans. As a result, an accurate estimate of the SI is the key in developing dependable forest management plans.

The site index depends on climate, geographic location, and tree species growing on the site [3]. Therefore, species-specific SI models are usually developed for even-aged pure stands. Morin et al. [4] simulated climate change impacts on the site productivity of European temperate forests across a large climatic gradient. They found that the impact varied across the gradient of current site conditions, irrespective of their composition. The forests on the warmest sites showed a decrease in productivity, while forest on sites with the coldest conditions experienced some productivity increases.

In mixed stands, however, two or more species growing at the same time, interspecies competition, and other factors may affect tree growth [5]. Consequently, SI models derived using the data collected from pure stands may not be suitable for mixed stands. Likewise, it may not be possible to develop SI models for each tree species growing in a mixed stand if the trees of targeted species are not appropriate (too small) to sample to develop a SI model for that tree species [6].

Where two species are growing in a stand and the SI of one species is known, it may be possible to determine the SI the other species in terms of the SI of the first species by developing a SI conversion equation. This equation represents a functional relationship between SIs of two species growing together in mixed stands. The relationship is usually linear and can be expressed using a linear equation. Stand dynamics of mixed species can also be modelled using this equation in a growth and yield application [7,8].

If there are two tree species growing together in mixed stands and both species have a few trees in the main canopy in each stand, the SI conversion equation between these species is developed by sampling dominant and codominant trees from both species [5,6,7,8]. If site indices of tree species A and B are represented by SI_A and SI_B, then their relationship can be expressed as follows:

SI_A = a + b SI_B + ε

(1)

where a and b are regression coefficients to be estimated, and ε is the error term.

Carmean et al. [9] synthesized findings from six studies that examined site index relationships across 24 hardwood and conifer species in the United States and Canada. They also derived site index conversion equations for trees grown in naturally regenerated mixed stands. Similarly, Nigh [5] utilized estimated site index values to formulate conversion equations for four conifer species in British Columbia, Canada.

There are mainly 3 types of naturally originated mixed stands in Ontario, Canada: (1) jack pine (Pinus banksiana Lamb.) and black spruce (Picea mariana Mill. B.S.P.), (2) black spruce and trembling aspen (Populus tremuloides Michx.), and (3) jack pine and trembling aspen. Climate effects on site productivity were analyzed, and climate sensitive stand height/site index models were developed for the first two types of these stands [10,11]. Site index conversion equations were also developed for the jack pine and black spruce mixed stands [10].

Climate-sensitive stand height/SI models and SI conversion equations are not available for naturally regenerated jack pine and trembling aspen mixed stands. Therefore, this study is focused on several objectives: (1) investigating the effects of climate on the site productivity of jack pine and trembling aspen trees in naturally regenerated mixed stands in Ontario, Canada; (2) developing models for stand height/SI for these species by integrating climate variables; (3) assessing the potential impact of future climate change scenarios on the stand height growth of both species; (4) comparing stand height growth between mixed and pure stands; and (5) establishing SI conversion equations for trembling aspen and jack pine grown in naturally originated mixed stands. These conversion equations are crucial for climate change adaptation and mitigation strategies as well as for wood supply modeling to estimate the SI of a new species in relation to the SI of currently growing species.

2. Materials and Methods

2.1. Height and Age Data

Stem analysis is a technique of examining the growth rings of sections from a tree trunk. It enables the past growth history of a tree stem to be reconstructed. Therefore, stem analysis was used to collect the height–age pair data used in this study. Thirty-one even-aged naturally regenerated mixed stands were sampled from across the study species’ ranges in Ontario (Figure 1). The longitude and latitude of these stands ranged from –79.7154° to –95.0769° and 47.5780° to 50.1767°, respectively. Similarly, the elevation covered from 212 to 483 m. These stands were dominated by trembling aspen (27.65%) and jack pine (23.37%). There was also a presence of other species in these stands. The proportions of other species were 16.46% black spruce, 14.35% white birch (Betula papyrifera), 13.5% balsam fir (Abies balsamea (L.) Mill.), and the remainder (4.67%) was a combination of white spruce (Picea glauca), white pine (Pinus strobus), red pine (Pinus resinosa), and red maple (Acer rubrum).

Jack pine and trembling aspen trees were sampled by establishing three 100 m² circular temporary sample plots at each site. One non-veteran largest diameter tree was sampled from each species from each plot for stem analysis. A veteran tree generally consists of a large trunk with one or more defects including large cavities, decay, and fungi on the trunk. Only trees without visible abnormalities, such as forks, broken or dead tops, and major stem injuries were selected for sampling. Stem analysis data were obtained by following the procedure described by Sharma [10]. The procedure described by Sharma [10] includes sampling trees, cutting disks from sampled trees, and conducting ring analysis. Therefore, readers are referred to those studies for further information on the procedure of collecting the height–age pair data used in this study.

The annual height growth for each tree of each species was estimated to gather height–age data (growth series). To ensure that the sampled trees had not experienced early height growth suppression, dieback, or breakage, the height of each tree was plotted against its age for each species. No defects were detected in any of the sampled trees. Consequently, data from all 186 trees (93 jack pine and 93 trembling aspen) were used to analyze and model the impact of climate on site productivity.

Generally, the height growth of trees before reaching breast height is inconsistent (erratic) in naturally originated mixed stands. Therefore, measured tree height from breast height and age from breast height (referred to as breast height age, BHA) were used in this study. Due to the erratic growth before breast height, site scale observations could not be made by combining growth series from the three plots at each site. Combining the series would average the height growth across years because the three trees from different plots at a site would reach BHA in different calendar years. Since climate variables are tied to specific calendar years, the climatic effects on height growth could not be analyzed if the series from the three trees across the plots were combined.

Summary statistics of the sampled trees for this study are presented in

Table 1. Summary statistics for site index (SI), total height, DBH, total age, breast height age (BHA), and climate variables for jack pine and trembling aspen trees used in this study. WQMT and WQTP represent the warmest quarter mean temperature and precipitation, respectively; MinTCP is the minimum temperature of coldest period; DQMT is the driest quarter mean temperature; and CMI_May and CMI_Dec are the climatic moisture indices of May and December, respectively. (N = number of samples and SD = standard deviation).

Variable	N	Mean	SD	Min	Max
Jack pine
Total age (year)	93	80.06	21.84	39.00	127.00
BHA (year)	93	76.53	21.63	36.00	124.00
Total height (m)	93	21.66	2.71	15.80	28.50
DBH (cm)	93	28.98	5.54	17.40	45.80
SI (m)	27	18.20	2.07	15.31	22.32
Trembling aspen
Total age (year)	93	77.84	19.77	41.00	129.00
BHA (year)	93	74.58	19.83	39.00	126.00
Total height (m)	93	22.03	2.71	16.64	31.50
DBH (cm)	93	29.88	5.46	19.00	43.50
SI (m)	29	18.53	1.83	14.96	21.97
Climate variables
WQMT (°C)	2335	16.17	1.33	12.40	19.90
WQTP (mm)	2335	247.33	46.67	115.00	415.00
MinTCP (°C)	2335	−24.81	3.18	−34.50	−16.60
DQMT (°C)	2335	−10.52	5.02	−20.50	14.90
CMIMay	2335	0.52	3.17	−7.01	10.24
CMIDec	2335	4.45	1.87	0.68	18.11

. The SI values in Table 1 indicate the average heights of three trees sampled from three plots at each site, measured at BHA 50 years for each species. The SI values are 27 for jack pine and 29 for trembling aspen. These values account for the fact that trees on four jack pine sites and two trembling aspen sites were younger than 50 years BHA.

2.2. Climate Data

The sample plots in this study were not near established weather stations, so Canadian climate models were used to estimate climate variables for each site [12]. These models were generated from continuous climate grids using ANUSPLINE based on corrected Canadian weather station data [13,14] provided by many stations in Ontario. Estimates of average yearly and seasonal values of these variables at each plot location were calculated for each year starting when the sampled tree reached breast height to 2015.

A total of 68 variables were calculated, including mean, minimum, and maximum air temperatures and total precipitation, estimated for each month of the year, for each quarter (consecutive three-month periods), and annually. In addition, climate data included estimates for start, end, and length of the growing season and the sum of growing degree days using a base temperature of 5 °C. The growing season was defined as the length of time between the day after 1 March when the mean daily temperature was ≥5 °C for five consecutive days and the day after 1 August when the minimum daily temperature was ≤−2 °C. The 68 variables also included three site-related variables (longitude, latitude, and elevation). In addition to the climate variables, climatic moisture index (CMI), obtained by subtracting monthly potential evapotranspiration (PET) from mean monthly precipitation (MMP) (see [15]), was estimated for each month for each year.

Estimates for all climate variables were provided by Dan McKenney (Canadian Forest Service, Personal Communication, 2023). Seasonal and (or) annual values of climate variables for a calendar year were used to examine and model the climate effects on the annual height growth of each tree for that calendar year. All these climate related variables are listed in

Table 2. Climate variables evaluated in analyzing climate effects in this study.

1.	Longitude
2.	Latitude
3.	Elevation
4.	Mean diurnal temperature range (MDTR) (mean (period max-min))
5.	Isothermality (MDTR/annual temperature range)
6.	Temperature coefficient of variation
7.	Max temperature of the warmest period
8.	Min temperature of the coldest period
9.	Temperature annual range
10.	Mean temperature of the wettest quarter
11.	Mean temperature of the driest quarter
12.	Mean temperature of the warmest quarter
13.	Mean temperature of the coldest quarter
14.	Annual precipitation
15.	Precipitation of the wettest period
16.	Precipitation of the driest period
17.	Precipitation coefficient of variation
18.	Precipitation of the wettest quarter
19.	Precipitation of the driest quarter
20.	Precipitation of the warmest quarter
21.	Precipitation of the coldest quarter
22.	Julian day number at the start of the growing season
23.	Julian day number at the end of the growing season
24.	Number of days of the growing season
25.	Total precipitation for period 1
26.	Total precipitation for period 3
27.	GDD above base temperature for period 3
28.	Annual mean temperature
29.	Annual minimum temperature
30.	Annual maximum temperature
31.	Mean temperature for period 3
32.	Temperature range for period 3
33.	January mean monthly minimum temperature
34.	February mean monthly minimum temperature
35.	March mean monthly minimum temperature
36.	April mean monthly minimum temperature
37.	May mean monthly minimum temperature
38.	June mean monthly minimum temperature
39.	July mean monthly minimum temperature
40.	August mean monthly minimum temperature
41.	September mean monthly minimum temperature
42.	October mean monthly minimum temperature
43.	November mean monthly minimum temperature
44.	December mean monthly minimum temperature
45.	January mean monthly maximum temperature
46.	February mean monthly maximum temperature
47.	March mean monthly maximum temperature
48.	April mean monthly maximum temperature
49.	May mean monthly maximum temperature
50.	June mean monthly maximum temperature
51.	July mean monthly maximum temperature
52.	August mean monthly maximum temperature
53.	September mean monthly maximum temperature
54.	October mean monthly maximum temperature
55.	November mean monthly maximum temperature
56.	December mean monthly maximum temperature
57.	January mean monthly precipitation
58.	February mean monthly precipitation
59.	March mean monthly precipitation
60.	April mean monthly precipitation
61.	May mean monthly precipitation
62.	June mean monthly precipitation
63.	July mean monthly precipitation
64.	August mean monthly precipitation
65.	September mean monthly precipitation
66.	October mean monthly precipitation
67.	November mean monthly precipitation
68.	December mean monthly precipitation
69.	Climatic moisture index

. Summary statistics of the climate variables that were significant in explaining the variability in stand height growth of trembling aspen and jack pine trees are displayed in Table 1.

2.3. Stand Height/Site Index Models

Tree height and age have a nonlinear relationship, so nonlinear mathematical expressions are typically used to describe it. These expressions are based on fractional and exponential functions. The Hossfeld IV and Chapman–Richards functions are the most commonly used functions for modeling the height–age relationship, representing fractional and exponential functions, respectively (see [16,17,18,19]). Various versions of these functions can be developed by applying different assumptions and initial conditions.

Sharma et al. [20] and Sharma and Parton [3,21,22] explored various forms (variants) of the Hossfeld IV and Chapman–Richards functions for four major commercial tree species in Ontario, Canada. They discovered that a particular variant of the Hossfeld IV function (Equation (2)) offered the best fit (lowest AIC [23]) and provided biologically plausible height estimates across productivity classes for these species. Consequently, this variant was used as the base function to analyze and model the climate effects on stand height growth in this study. The specific variant (model form) they used is provided below:

$$H_2={\displaystyle \frac{{\alpha }_0}{1-\left(1-\frac{{\alpha }_0}{H_1}\right){\left(\frac{A_1}{A_2}\right)}^{{\alpha }_1}}}+\epsilon$$

(2)

where H₁ and H₂ are stand heights (above breast height) at BHAs A₁ and A₂, respectively, α₀ and α₁ are parameters and ε is the error term. This model form is also known as McDill–Amateis growth function (see [24]). In general, α₀ represents the asymptote of the growth curve, and α₁ determines the shape. Here, α₁ is also known as the rate parameter. As described in previous studies, climate effects on stand height growth can be analyzed and modelled by expressing the parameters (α₀ and α₁) in terms of climate variables. Climate variables can also be incorporated in Equation (1) to develop SI conversion equations.

2.4. Model Fitting and Evaluation

Height–age measurements used in this study are organized hierarchically, with measurements nested within individual trees and trees nested within specific site. As a result, observations among sites are independent but measurements within a tree are dependent and correlated. The autocorrelation resulted by this data structure was addressed by using a mixed-effects modeling approach. Similarly, heteroscedasticity (if present) was addressed by specifying a variance function [25]. Two variance functions (power and exponential) were taken into consideration, and the function resulting in the best fit (smallest AIC value) was used in the model.

Initially, no climate variables were used to fit Equation (2) to the data from both jack pine and trembling aspen trees in SAS [26]. Climate variables were then incorporated into the model by categorizing them into two groups (temperature and precipitation) and adding one variable from each group at a time. From each group, the variable that was significant in the model (α = 0.05) and resulted in the lowest fit statistics (Akaike’s information criterion (AIC [23]) and MSE) was selected. Next, a climate variable from the second group was introduced into the model, one at a time, in the presence of the variable selected from the first group. The variable that resulted in the best fit statistics and was significant in the model was then selected from the second group. This process continued until no additional variable from either group was significant in the model for both species.

Site-specific variables were then incorporated into the model if they were significant and improved fit statistics in the presence of climate variables. Random effects were incrementally added to fixed-effects parameters if significant. To assess heteroscedasticity in the data, residuals (observed–predicted) were calculated for all one-year growth periods of each tree species and plotted against predicted heights.

Climate effects on future stand height growth were assessed by predicting stand heights for both species across four regions (represented by triangles) in northern Ontario, including eastern (near New Liskeard), central (near Geraldton), western (near Thunder Bay), and far western (near Manitoba border), as illustrated in Figure 1. These predictions were generated using a model incorporating climate variables under three representative concentration pathways (RCPs), which simulate different levels of warming (2.6, 4.5, and 8.5 Watts/m²) by the end of the century according to the Canadian model [12]. The projected climate variables deemed significant for defining parameters in the height growth model were used to evaluate climate impacts. Growth curves for the 79-year period from 2022 to 2100 were then generated for assessments.

Two separate SI conversion equations were necessary to estimate the site index (SI) of jack pine based on trembling aspen, and vice versa, because the climate variables influencing the SI varied between species. Therefore, distinct SI conversion equations were formulated for mixed stands of jack pine and trembling aspen. These equations incorporated climate variables that significantly affected the SI of each species. For this, the SI of trembling aspen was regressed against the jack pine SI and relevant climate variables, and vice versa. The SI values used in deriving SI conversion equations were determined using stand height at a breast height age (BHA) of 50 years. Climate variables used in these equations represented the average values over the 50-year growth period.

3. Results

3.1. Climate Effects on Stand Height Growth

Equation (2) was initially analyzed using height–age data from jack pine and trembling aspen trees. Random effects related to the rate parameter and asymptote, along with a variance function, were included in the model fitting process. Both random effects were found to be significant at the stand scale for both species. Both variance functions (power and exponential) showed high significance in the regression analysis. However, the power function applied to height resulted in a smaller AIC for both species. Therefore, the height growth model incorporating random effects and the power variance function is formulated as follows:

$$H_{ijk}={\displaystyle \frac{({\alpha }_0+b_{0i})}{1-\left(1-\frac{\left({\alpha }_0+b_{0i}\right)}{H_{ijl(k\ne l)}}\right){\left(\frac{A_{ijl(k\ne l)}}{A_{ijk}}\right)}^{{\alpha }_2+b_{1i}}}}+{\epsilon }_{ijk}$$

(3)

where H_ijk is the stand height at age A_ijk (k^th observations of tree j at stand (site) i), H_ijl is the stand height of the same tree at age A_ijl at the same site (l^th observations of tree j at site i and k ≠ l), and b_0i and b_1i are stand-scale random effects. These random effects are connected to α₀ and α₁, respectively, and are normally distributed with mean zero and variances σ₀² and σ₁², respectively (i.e., b_0i ~ N(0, σ₀²) and b_1i ~ N(0, σ₁²)). In addition, these are independent of ε_ijk ~ N(0, σ²H_ijk^δ). Delta (δ) is the power to the height (H_ijk) of the variance function. Other variables have been defined previously.

Introducing random effects in Equation (2) improved fit statistics (AIC and MSE) for both species. However, the estimated fixed effects parameters were not biologically consistent for jack pine; specifically, both the asymptote and rate parameters were overestimated when random effects were included. In contrast, for trembling aspen, the estimated values of both fixed effects parameters in the presence of random effects were reasonable and biologically consistent. Therefore, Equation (3) without random effects was employed to analyze and incorporate climate effects on stand height growth for jack pine trees. The estimated parameters for both species are presented in

Table 3. Parameter estimates, standard error (SE), MSE (σ_e²), variance of b₀ (σ₀²), variance of b₁ (σ₁²), covariance of b₀ and b₁ (σ₀ σ₁), power variance (δ) to the height, and Akaike’s information criterion (AIC) for Equation (3) for jack pine and trembling aspen trees grown in natural mixed stands in Ontario, Canada.

Parameters	Jack Pine		Trembling Aspen
	Estimates	SE	Estimates	SE
α0	34.8151	0.3006	36.7482	0.8591
α1	1.1605	0.0083	1.1526	0.0308
σe2	0.2751	--	0.0727	--
σ02	--	--	120.0100	--
σ12	--	--	0.0251	--
σ0σ1	--	--	−1.0689	--
δ	−1.4094	0.0351	−0.6840	0.0330
AIC	−13,249	--	−9437	--

.

Equation (3) initially lacks climate variables. To explore the influence of climate on stand height growth, the parameters (α₀ and α₁) in Equation (3) were reformulated in terms of climate variables. As outlined earlier, climate variables were categorized into two groups: temperature and precipitation. Site-specific variables such as elevation, longitude, and latitude constituted a third group (site). Each fixed effects parameter (α₀ and α₁) in Equation (3) was subsequently modeled as a function of all 69 variables from each group individually.

When analyzed individually, many temperature- and precipitation-related variables showed significance in explaining variations in both parameters. However, the asymptote (α₀), when expressed in terms of warmest quarter total precipitation (WQTP), yielded the best fit (lowest AIC and MSE values) for both species. In the presence of WQTP (related to precipitation), the rate parameter, expressed as a function of warmest quarter mean temperature (WQMT), provided the best fit for both species. With both WQTP and WQMT included, however, none of the other climate or site-related variables were significant in the model for either species. Random effects associated with the asymptote were not significant in the presence of climate variables for trembling aspen. Therefore, the final height growth models incorporating climate variables are formulated as follows:

Jack pine

$$H_{ijk}={\displaystyle \frac{{\alpha }_0+{\alpha }_2{WQTP}_i}{1-\left(1-\frac{{\alpha }_0}{H_{ijl(k\ne l)}}\right){\left(\frac{A_{ijl(k\ne l)}}{A_{ijk}}\right)}^{\left({\alpha }_1+{\alpha }_{3 }{WQMT}_i\right)}}}+{\epsilon }_{ijk}$$

(4)

Trembling aspen

$$H_{ijk}={\displaystyle \frac{{\alpha }_0+{\alpha }_2{WQTP}_i}{1-\left(1-\frac{{\alpha }_0}{H_{ijl(k\ne l)}}\right){\left(\frac{A_{ijl(k\ne l)}}{A_{ijk}}\right)}^{\left({\alpha }_1+b_{1i}+{\alpha }_{3 }{WQMT}_i\right)}}}+{\epsilon }_{ijk}$$

(5)

where WQTP_i and WQMT_i are the warmest quarter total precipitation and mean temperature, respectively, at site i; α₂ and α₃ are coefficients of these climate variables; and all other variables are as previously defined. Estimated values of parameters for Equations (4) and (5) are listed in

Table 4. Parameter estimates, standard error (SE), MSE (σ_e²), variance of b₁ (σ₁²), power variance (δ), and Akaike’s information criterion (AIC) for Equation (4) for jack pine and (5) for trembling aspen trees grown in naturally originated mixed stands in Ontario, Canada.

Parameters	Jack Pine		Trembling Aspen
	Estimates	SE	Estimates	SE
α0	39.6855	0.7406	41.6294	1.3594
α1	1.9344	0.0563	1.74140	0.0964
α2	−0.02135	0.0026	−0.02202	0.0050
α3	−0.04903	0.0035	−0.01580	0.0057
σe2	0.2605	--	0.08811	--
σ12	--	--	0.1661	--
δ	−1.4057	0.0310	−0.7356	0.0318
AIC	−13,553	--	−9090	--

. Since introducing climate variables significantly improved fit statistics for both tree species, Equations (4) and (5) can be used to describe climate effects on the stand height growth of jack pine and trembling aspen, respectively, grown in naturally originated mixed stands.

Coefficients of WQTP and WQMT were both negative for both species. As a result, increasing WQTP will decrease the asymptote for both species. Similarly, increasing WQMT will negatively affect the rate of height growth of both jack pine and trembling aspen. Once the model was fit with the variance function, no further heteroscedasticity was detected for both species. This was confirmed by plotting residuals against the predicted values of stand heights. Estimates for the power of height (weight) were negative. The AIC values were decreased significantly by incorporating weight into the model for both tree species.

Climate effects on future stand height growth was assessed by predicting stand heights using Equations (4) and (5) for jack pine and trembling aspen, respectively. These predictions were made across four sites in Ontario under three emissions scenarios (RCPs) for the growth period from 2022 to 2100. Initial values of heights used here were derived from average heights of species-specific sampled trees across all 31 sites at a breast height age (BHA) of 5 years. Additionally, yearly projections of climate variables (warmest quarter total precipitation, WQTP, and warmest quarter mean temperature, WQMT) under three climate change scenarios were utilized in estimating these heights. The height-age pairs were then used to generate height–age curves (Figure 2 and Figure 3).

Climate effects were apparent on stand height growth from the beginning of the growth period for both species under all emission scenarios at all locations. The effects were negative for jack pine but positive for trembling aspen. For jack pine, the negative effects under RCP 2.6 were almost identical to those under RCP 4.5 at all locations. Height growth curves under RCPs 8.5 and 2.6 separated from each other approximately after age 50 years in the western, central, and far western locations. In the east, however, these curves separated much earlier (around 25 years), and the separation at this site was more pronounced than those at other locations. The most and least pronounced climate effects on height growth for jack pine were in the far west and central locations, respectively. Stand heights at BHA 81 years under RCP 8.5 scenario at central, western, eastern, and far western areas were shorter by 11.15, 16.66, 19.66 and 23.33%, respectively, than those estimated using Equation (3) that does not include any climate variables (Figure 2).

In the case of trembling aspen stand height growth, the overall effects of climate were positive under all emission scenarios at all locations (Figure 3). At a particular location, height growth curves under all three scenarios were almost indistinguishable from one another. The variation in climate effects on the height growth was not pronounced across locations. Stand heights at a BHA of 81 years under RCP 8.5 at western, far western, eastern, and central areas were taller by 12.22, 12.44, 12.78 and 14.37%, respectively, than those estimated using Equation (3) without climate variables.

Finally, the growth of stand heights in jack pine and trembling aspen trees, cultivated in both natural mixed and pure species stands, was compared. This involved estimating stand heights across different stand ages using the model without climate variables (Equation (3)) and models from previous studies with the same form. For jack pine, heights were estimated using Equation (3) for trees grown alongside trembling aspen (this study), Equation (3) as published by Sharma and Reid [27] for naturally occurring pure stands, and Equation (8) as published by Sharma [10] for mixed stands alongside black spruce trees.

For trembling aspen, these heights were estimated using Equation (3) for mixed stands grown with jack pine (this study) and Equation (2) published by Sharma [11] for mixed stands grown with black spruce trees. A stand height growth model for trembling aspen grown in pure natural stands was not available. Estimated heights were plotted against stand age for both tree species for all stand types (Figure 4). Average heights at a BHA of 5 years were used as initial heights in generating height growth curves for both species.

For jack pine, there was no difference in heights between pure and black spruce mixed stands until around BHA 30 years. The difference in height of jack pine mixed with trembling aspen, however, was obvious after an approximate BHA of 10 years. Thereafter, height in trembling aspen mixed stands outpaced that in pure and mixed with black spruce stands, with the difference increasing over time. For trembling aspen, height in stands mixed with black spruce was consistently lower than that in stands mixed with jack pine trees from the beginning. This indicates that interspecies competition is different for different tree species as found in [4].

3.2. Site Index Conversion Equations

The site index (height measured at an index age of 50 BHA) of trembling aspen trees was regressed linearly against the SI of jack pine trees grown in the same plot, and vice versa. The regression model is expressed as follows:

y_ij = β₀ + β₁ x_ij + ε_ij

(6)

where y_ij is the SI of trembling aspen (jack pine) at site i and plot j, β₀ and β₁ are the parameters to be estimated, and x_ij is the SI of jack pine (trembling aspen) tree at site i and plot j. In this case, individual plot SI values were used to develop SI conversion equations so that climate variables could be included in the model. As mentioned earlier, trees from three plots at a site did not reach breast height during the same calendar year.

Equation (6) was fit without climate variables first for both species. Estimated parameters are displayed in

Table 5. Estimated parameters (standard errors in parentheses) and fit statistics (MSE (σ_e²), Akaike’s information criterion (AIC), and R²) for Equations (6)–(8) for jack pine and trembling aspen trees from Ontario, Canada.

Parameters	Jack Pine		Trembling Aspen
	Equation (6)	Equation (7)	Equation (6)	Equation (8)
β0	6.1024 (2.0068)	4.14456 (1.8105)	9.3209 (1.523)	35.9208 (5.3729)
β1	0.6529 (0.1078)	0.35609 (0.0172)	0.5031 (0.083)	0.2062 (0.0912)
β2	--	−1.67383 (0.3781)	--	1.4858 (0.2649)
β3	--	1.88047 (0.2996)	--	−1.0911 (0.2102)
β4	--	--	--	1.1423 (0.2054)
σe2	3.9434	2.6283	3.0384	1.9191
R2	0.3285	0.5644	0.3285	0.5928

. In this case, the coefficient of determination (R²) for both species was 32.85% (Table 5). Climate variables were then included to the right-hand side of Equation (6) and fit to SI data from both species. The following models were selected based on R², AIC, and VIF (variance inflation factor).

Jack pine

PJSI_ij = β₀ + β₁ PTSI_ij + β₂ (CMI_May)_ij + β₃ (CMI_Dec)_ij + ε_ij

(7)

Trembling aspen

PTSI_ij = β₀ + β₁PJSI_ij + β₂ MinTCP_ij + β₃ DQMT_ij + β₄ (CMI_Dec)_ij + ε_ij

(8)

where PJSI_ij and PTSI_ij are the SIs of jack pine and trembling aspen, respectively; (CMI_May)_ij is the May climatic moisture index; (CMI_Dec)_ij is the December climatic moisture index; MinTCP_ij is the coldest period minimum temperature; DQMT_ij is the driest quarter mean temperature at site i and plot j; and β₁–β₄ are regression coefficients. Other variables are as previously defined.

Equations (7) and (8) were fitted using PROC REG in SAS. The introduction of climate variables significantly improved both fit statistics (R² and MSE) (Table 5). All climate variables were highly significant (p < 0.01) in both equations. The coefficient of determination (R²) increased from 0.3285 to 0.5644 for jack pine and from 0.3285 to 0.5928 for trembling aspen trees by including climate variables.

The introduction of climate variables into the model also improved the distribution of residuals (observed—predicted) (Figure 5). The residual plot for jack pine was very similar to that for trembling aspen. The correlation between climate variables used here was negligible as the maximum value of VIF was less than 3.33 for CMI_Dec in Equation (7) and 3.47 for MinTCP in Equation (8), and the rest were less than 3.0 for both models.

4. Discussion

The productivity of a site for growing trees is affected by environmental conditions including climate [28]. Changes in temperature and precipitation patterns have notably impacted the site productivity of various tree species, such as red pine, jack pine, white pine, black spruce, and white spruce, whether grown in plantations or naturally occurring pure and mixed stands across Ontario [3,10,11,20,21,22]. The direction (positive/negative) of these effects varied depending on the tree species and geographical location. Furthermore, the extent of climate influence differed among tree species and exhibited variations from east to west and north to south, even within the same stand type (plantations versus natural stands) and between pure and mixed stands. Additionally, climate variables affecting productivity also showed variability across different species and stand types [3,10,11,20,21,22,27].

In monospecific plantations, variables associated with precipitation and temperature—specifically, total precipitation of growing season (TPGS) and mean temperature of the growing season (MTGS)—played a significant role in influencing stand height growth for jack pine and black spruce trees [20]. Conversely, for planted red pine, only a temperature-related variable (MTGS) was found to be significant in the stand height growth model [22]. Similarly, the mean temperature of warmest quarters (MTWQ) and total precipitation of warmest quarters (TPWQ) were influential in explaining stand height growth variability for planted white spruce [21]. In contrast, average diurnal temperature range (ADTR) was more significant in describing the variation in stand height growth for white pine than other climate variables [3].

In naturally regenerated mixed stands, only one temperature-related variable (MTGS) was sufficient to explain the variability in stand height growth for both black spruce and jack pine trees [10]. Similarly, the average diurnal temperature range (ADTR) was the predominant climate variable influencing height growth models over others for both trembling aspen and black spruce trees [11]. However, in this study, both temperature-related variables (MTWQ) and precipitation-related variables (TPWQ) explained the most significant portion of variation in stand height growth for both jack pine and trembling aspen trees compared to other climate factors. These findings highlighted that climate effects on site productivity varied with stand type.

The results in this study showed that the interaction between jack pine and black spruce was negative (competition) but was positive (facilitation) between jack pine and trembling aspen as previously described in [29,30]. As far as climate effects are concerned, increasing the temperature of the warmest quarter could decrease the productivity of both species. Morin et al. [4] also reported that the site productivity of temperate European forests was strongly impacted by climate change, and the impact varied in the direction as well as in magnitude. The forests on the warmest sites showed a decrease in productivity, while forest on sites with the coldest conditions resulted in some increase in productivity. Their study also showed that indirect effects related to species composition on site productivity was also important in addition to direct and indirect effects of climate change. This corroborates the results found in this study.

Sharma [10] utilized observed site index (SI) values to establish SI conversion equations for jack pine and black spruce in naturally regenerated mixed stands. This was achieved by regressing the SI of one species against the SI of the other species while incorporating climate variables. The significant climate variables differed between species, except for the mean monthly minimum temperature of December, which remained consistent. Additionally, Sharma [11] analyzed SI data from measured values of trembling aspen and black spruce trees in naturally occurring mixed stands. He found no clear relationship between the SIs of trembling aspen and black spruce trees; hence, no SI conversion equations were reported.

In this study, the correlations between SIs of jack pine and trembling aspen in naturally regenerated mixed stands was examined using SI values derived from stem analysis. The SI of jack pine was found to be linearly associated with the SI of trembling aspen, and vice versa. Incorporating climate variables into the SI conversion equations strengthened this linear relationship between the two species. The climatic moisture index for December (CMI_Dec) was the only climate variable significant in both conversion equations. Additionally, CMI_Dec and the climatic moisture index for May were significant in the SI equation for jack pine. In contrast, for trembling aspen SI, in addition to CMI_Dec, the minimum temperature during the coldest period and the mean temperature of the driest quarter were significant.

The SI conversion equations that Nigh [5] developed for mixed stands were obtained by fitting these equations to SI data estimated using models fit to data gathered from naturally originated pure stands. However, findings from the aforementioned studies indicated that the growth patterns in mixed stands often diverge from those observed in pure stands. This contrast was particularly noticeable for black spruce and jack pine cultivated in both mixed and pure stands [10,11]. Consequently, SI models formulated for mixed stands using data derived from pure stands might lead to inaccurate conclusions.

In this study, stand heights were projected using the initial values at a breast height age (BHA) of 5 years to assess climate impacts for both species (see Figure 2 and Figure 3). Stand heights estimated from initial values at BHAs less than 5 years did not align with measured heights (significantly taller than observed heights in those stands). Hence, for future projections, it is recommended to base height estimations on initial heights measured at least at 5 years from BHA for both jack pine and trembling aspen trees cultivated in naturally occurring mixed stands.

5. Conclusions

Conversion equations for site index (SI) were established for jack pine and trembling aspen trees grown in naturally regenerated mixed stands in Ontario, Canada. Initially, the influence of climate on stand height growth for both species was examined, subsequently developing stand height growth models that integrated climate variables. The McDill–Amateis growth function served as the analytical tool to evaluate climate impacts and construct climate-sensitive models for the stand height growth of these tree species in naturally regenerated mixed stands.

Without climate variables, the correlation between SI values of jack pine and trembling aspen was modest (R² = 0.3285). However, incorporating climate variables into the models substantially enhanced this relationship between the SI values of these tree species (R² = 0.5644 for jack pine and 0.5928 for trembling aspen).

The height growth models developed in this study allow for the estimation of stand heights and site indices (SIs) for jack pine and trembling aspen trees grown in naturally regenerated mixed stands. Additionally, SI conversion equations enable the estimation of one species’ SI based on the SI of another species and environmental factors. If data on climate variables are unavailable, models that do not incorporate these variables can still be applied to estimate the stand height and site index of both species. Since site productivity is affected by climate and SI is the main driver of most of forest growth and yield models, climate-sensitive SI models should be developed for all major commercial tree species.