Effects of site productivity on individual tree maximum basal area growth rates of Eucalyptus pilularis in subtropical Australia

Inventory data were available from 96 plots of even-aged, monoculture, tall-open forests of Eucalyptus pilularis Smith, aged 2–63 years, growing in sub-tropical regions along the east coast of Australia. A model was developed relating the maximum possible stem basal area growth rate of individual trees to their stem basal area. For any tree size, this maximum increased as site productivity increased. However, the size at which this maximum occurred decreased as productivity increased. Much research has shown that, at any stand age, trees of a particular stem basal area are taller on more productive sites than on less productive ones. Taller trees incur greater respiratory costs to ensure maintenance of the photosynthetic capacity of their canopies; this reduces their growth rates. It was concluded that trees with larger basal areas will have the maximum possible growth rate on a less productive site, whilst trees with smaller basal areas will have the maximum possible on a more productive site. The model developed may constitute the first stage of a complete individual tree growth model system to predict wood yields from these forests.


Introduction
One approach to the development of models to predict growth rates of individual trees in forests is to first construct a model predicting the maximum possible growth rate (sometimes termed potential or optimum growth rate) that a tree might have in relation to its size at the time growth is occurring (Smith et al. 1992;Pretzsch and Biber 2010;Weiskittel et al. 2011;Pommerening et al. 2022).Terms are then added to the model to predict the reductions from this maximum that are a consequence of both the particular environmental circumstances of the site on which the tree is growing and the competition it faces from its neighbours for the resources it requires for growth (light, water and nutrients).
The present work used growth data obtained from inventory plots of Eucalyptus pilularis Smith growing in sub-tropical regions along the east coast of Australia.Native forests of this species occur as tall-open forests that are generally even-aged and are often monocultures.The species has been of major commercial importance for timber production for many years and is sometimes grown in plantations.West (2021) developed a model relating the overall maximum stem basal area growth rates of individual trees to their basal areas.An objective of the present study was to extend that model to determine how those growth rate maxima were affected by changes in the productivity of the site on which the forest was growing.A second objective was to then relate the effects of both tree size and site productivity on growth rates to known physiological behaviour.
The model developed here should prove a useful starting point for an individual tree basal area growth rate model to which the effects of inter-tree competition and probability of death can be added subsequently to render a complete model useful for forest management.

Data Sources
The data used here had been collected by the New South Wales and Queensland government forest services during 1931-1989 from 96 inventory plots located between 26-29° S and 152-153° E; this is only part of the total area of natural occurrence of E. pilularis (Boland et al. 2006).It should be appreciated that these data do not represent a formal sample of these forests but were those data available immediately over the period they were being collated (1983)(1984)(1985)(1986)(1987)(1988)(1989)(1990)(1991) by Mattay and West (1994).
The plots were selected to be in even-aged stands and were essentially monospecific, with only rare inclusions of individuals of other species.Forest ages ranged over 2-63 years at the time of measurement.Plots varied in area over 0.04-0.2ha and were measured 3-28 times at intervals of 0.4-5 years.At each measurement of a plot, diameters at breast height (1.3 m) over bark of live trees were recorded and, often, the heights of at least a sample of trees.In the entire data set, the largest diameter was 80 cm and the tallest individual was 55 m.
For every tree that survived between plot measurements, its basal area growth rate was determined as its change in basal area divided by the length of the growth period.Where data have come from growth periods longer than several years, others have determined values of basal area growth rate that incorporate information about the trend of growth rate with age over the growth period concerned (McDill and Amateis 1993;Cao 2000;Cao et al. 2002;Yue et al. 2012).
In the present study, 90% of the growth rates were determined from growth periods of three years or less; it was felt there would be little advantage from using these more complex methods here.A total of 1,042 growth periods were measured over the 96 plots, giving 35,876 individual tree basal area growth rates.
Where heights had been measured, stand dominant height (average height of the 50 largest diameter trees per ha in a plot) was determined.Based on the dominant height at an age closest to 20 years, site productivity of each plot was determined as a site index, being the predicted dominant height at an index age of 20 years; site indices were estimated using functions developed by West and Mattay (1993).Site indices in the data set ranged over 13-40 m.

Selecting growth rate maxima
Models that aim to predict maximum possible individual tree growth rates on a site assume that data are available for trees growing free of competition.Only if the trees are open grown, is it possible to be certain of this (Smith et al. 1992).Given broad-scale plot data such as available here, it is impossible to identify exactly which trees conform to this requirement.However, when data are available for thousands of growth rates, perhaps it may be assumed that there are at least a few individuals growing essentially under those circumstances.Most models of this type have been developed assuming this is the case.
Thus, the first step in model development was to reduce the available data to leave only those for trees which had high growth rates for their size and which it might be hoped were free growing.This has been done in other studies with different species around the world (Coomes and Allen 2007;Pretzsch and Biber 2010;Pommerening et al. 2011;Bošelǎ et al. 2013;Pommerening and Särkkä 2013;Pommerening and Maleki 2014;Häbel et al. 2019;Tian et al. 2021) using quantile regression (Koenker and Hallock 2001;Cade and Noon 2003).Other methods are available to achieve equivalent results (Bi and Turvey 1997;Bi et al. 2000;Zhang et al. 2005;Tian et al. 2021).For this study, the method of Bi and Turvey (1997) was used, termed here the 'class maxima' method.Note that, incorrectly, it was termed quantile regression in West (2021West ( , 2022)).This method was considered appropriate for the present work because it uses data that include, unequivocally, the maximum growth rates observed for the tree sizes encountered in the data set.
To use the class maxima method, the data set of 35,876 growth rate observations was split into 10 subsets, each containing trees from plots within a particular site index range.In the entire data set, trees were spread somewhat unevenly in numbers across plots with different site indices.However, attempts were made to choose site index classes that were roughly equal in width and contained substantial numbers of individual tree growth rates.Table 1 shows these site index classes and the number of observations in each data subset.Each subset was then split into 50 equal sized basal area classes and the observation with the highest growth rate in each size class was selected.Occasionally there were size classes which contained no observations and were ignored.A total of 460 observations were selected by this process over the 10 data subsets.
Figure 1 shows a scatter plot of the stand stocking densities against stand basal areas for the 1,042 plot measurement periods in the full data set.Superimposed on that is a corresponding scatter plot for the 460 observations that made up the 10 class maxima data subsets of maximum tree growth rates; many of those observations came from the same growth periods of the same plots and so share the same data point.The results suggest that the class maxima data tended to come from stands with lower stand basal areas and stocking densities, circumstances that are the most likely to include data for trees that are suffering little or no reductions in growth rates due to competition, as discussed previously.Note also, that quantification of the effects of inter-tree competition in these forests (West 2022) suggested there are indeed instances where trees were growing free, or virtually so, of competition from their neighbours.
Scatter plots in Fig. 2 show the basal area growth rates against basal areas of individual tree data selected by this process for four of these class maxima data subsets; these four were chosen as examples to reflect the overall range of plot site indices in the whole data set.It was assumed that these growth rates in each subset were from trees growing with little competition from other trees and thus reflected the maximum growth rate possible for a tree of its size on a site with a productivity within the site index range used for each class.The supplementary material (Fig. S1) with this paper includes these scatter plots for all 10 class maxima data subsets.

Modelling Model development
Model development was a two-stage process.The first stage involved the selection of a model function that would describe the selected class maxima data in any one of the 10 data subsets.The well-known Chapman-Richards function seemed suitable for this purpose (Pienaar and Turnbull 1973;Smith et al. 1992;Zeide 1993;West 2021) (Eq.1): where ΔB (m 2 a −1 ) is the basal area growth rate of an individual tree with basal area B (m 2 ), at the start of the growth period concerned, and α, β and γ are parameters.The nonlinear Model (1) was fitted to the data of each of the 10 data subsets by least-squares regression using the NLIN procedure of the SAS statistical package. 1 The parameter estimates for each data subset are shown in Table 1; the fit to the data of this model for each example data set is shown in Fig. 2.
A second-stage model was then developed to attempt to determine how site productivity, as estimated by stand site index, affected tree growth rates.Scatter plots were drawn of the first-stage parameter estimates against the average site index (S, m) of the plots from which the data had been   included in each data subset (Table 1).Examination of these scatter plots suggested that, of the three parameters, the value of the parameter β was the one most closely related to stand site index.The scatter plot of those data is shown in Fig. 3a.It appeared that the model could be used to describe the relationship, albeit rather poorly; the ordinary least-squares straight line fit of the (2a) = exp 1 + 2 ln S  West (2022) data, using the logarithmic expression of Model (2a), had an r 2 value of only 21%.It was then recognised that both parameters α and γ were related closely to β over the 10 data subsets; the corresponding scatter plots are shown in Figs.3b and 3c.It then appeared that appropriate secondstage models for those two parameters might be: and Equation 2a was substituted into each of Eqs.2b and 2c and the results then substituted into Eq. 1.This yielded a rather complex, nonlinear model involving the six parameters α 1 , α 2 , β 1 , β 2 , γ 1 and γ 2 , where each of β 1 and β 2 appeared three times at different positions in the model.This model may be written more fully as: where The NLIN procedure of the SAS package readily fitted this complex model to the 460 observations pooled over the 10 class maxima data subsets.The parameter estimates with their standard errors are shown in Table 2.
A potential problem with fitting a model to this pooled data set was that, whilst there may be random variation in the observed data of any one data subset for any one site index class, its average behaviour might differ from that of other data subsets of the same class that might have occurred, by chance, in the data.Under these circumstances, the residuals from the fitted model in any one data subset may no longer be independent of each other, a principal requirement of least-squares regression analysis.This correlation (covariance) between the residuals within a data subset may not lead to bias in the parameter estimates obtained for the model, but it must be accounted for if their standard errors are to be estimated without bias.This is the same problem as encountered in data sets that contain multiple measurements within individual sampling units (Cole and Grizzle 1966;Dielman 1983;West et al. 1984;Fitzmaurice et al. 2011).As West and Ratkowsky (2022) noted, various solutions to this problem have been suggested but not all are applicable with nonlinear models or can be relied upon on under all circumstances.
In the present case, Model (3) appeared to fit the data reasonably well.The fit in each of the four example data subsets is shown as the solid lines in Fig. 2. Inevitably, since such fitted lines represent predictions made using the average site index of the plots that made up the data in each subset, there appear to be some deviations from the fit to the data for each of the example data sets (the solid lines).However, scatter plots of the residuals against the fitted values from Model (3) did not show any marked pattern nor any tendency of a positive or negative bias in any subset.This suggested that the estimate of the covariance matrix of the parameter estimates was unlikely to have any substantial bias.The corresponding 95% confidence limits about estimates of basal area growth rates using the pooled model are shown as the more complex dashed f , S = exp 1 + 2 ln S lines in each of the examples in Fig. 2.These confidence limits were determined using the method of Box (1971), as shown in Eq. 9.6 of Ratkowsky (1983), and the estimate of the covariance matrix of the parameter estimates provided by the NLIN procedure of the SAS package.The fit to the data of Model (3) for all 10 class maxima data subsets is shown in the Supplementary material.
As an alternative to this use of nonlinear regression to fit Model (3), it is perhaps inviting to fit the set of Models (2a-c) jointly using the 3SLS method of seemingly unrelated regression (SUR) (Zellner 1962(Zellner , 1963;;Srivastava and Giles 1987).It was found that whilst this method gave an adequate fit to the data, it grossly overestimated the confidence limits about the fitted lines.This was probably a consequence of having only 10 observations to fit the model with SUR, whereas the nonlinear regression (Model 3) was based on 460 observations.An earlier model relating maximum tree basal area growth rates to tree basal areas and site productivity was developed by West (2022) using largely the same inventory data.The fit to the data of that model is shown in Fig. 2 with long dashed lines.It fitted the data well for stands with mid-level site indices (Fig. 2b and c), but underestimated growth rates in stands with higher site indices (Fig. 2a) and over-estimated in stands with lower site indices (Fig. 2d).That model was less complex in form than the present model and this probably led to its somewhat poorer fit to the data.

Application of model
Figure 4 shows predictions from Model (3) of individual tree maximum basal area growth rates in relation to basal areas for each of a set of five site indices.It shows clearly how maximum growth rate for trees of any basal area increased with increasing site index, consistent with the trends in the data for the four site index classes in Fig. 2.
Within any one stand, the model predicts that maximum growth rate increases progressively with tree size, to an overall maximum for the stand and then declines progressively as tree size increases further.This is the form predicted by the Chapman-Richards function (Model 1).With straightforward calculus, it is easy to show that the basal area (B m , m 2 ), at which that overall maximum occurs in any stand, may be predicted as: For a stand of any site index S̄ (m), Models (2a-c) may be used to give values for α, β and γ for use in Eq. 4.This was used in drawing the dashed line in Fig. 4. The corresponding overall maximum growth rate for that site index may then be determined by applying B m for B in Model (1).( 4) The dashed line (Fig. 4) shows how B m tends to increase progressively as stand site productivity declines.Over the site index range 20-40 m, the progressive increase in B m corresponds to a progressive increase in stem diameters over the range 22-34 cm.The corresponding overall maximum growth rate declines progressively as stand site productivity declines.

Discussion and conclusions
As discussed in earlier work (West 2021(West , 2022)), the Chapman-Richards function as displayed here is consistent with the well-known phenomenon that stand growth rates of evenaged forests increase during their early years, but within a few decades, reach a maximum and decline progressively thereafter.Various reviews and texts describe this phenomenon (Gower et al. 1996;Ryan et al. 1997;Davis et al. 2001;West 2015;Binkley 2023).These trends are evident in whatever measures of tree size are used, measures such as biomass, wood volume or characteristics correlated with these, such as stem basal area.During early stand development, individual tree growth rates increase with size, hence age.Stand growth is the sum of individual tree growths, so as trees become larger, stand growth rate will become correspondingly bigger.However, after some years, some trees will reach a size at which their growth rates decline, consistent with the trends in maximum growth rates shown in Figs. 2 and 4.Over time, more and more trees will reach this size and stand growth rates will continue to decline, the decline increasing steadily with increasing tree size, that is, with increasing stand age.
Since the mid-twentieth century, forest research has sought to explain this growth behaviour.Consideration of various alternatives led Binkley (2023) to conclude that the most likely cause was a reduction with age (hence size) of the 'efficiency' of tree growth, i.e., the amount of growth per unit of light, water and nutrients that a tree obtains from the environment.A theory consistent with this, which has received considerable support, is the 'hydraulic limitation theory'.This proposed that water stress in leaves increases as taller trees must transport water from roots to leaves to greater heights; this was assumed to reduce the photosynthetic capability of their canopies (Ryan et al. 1997(Ryan et al. , 2006;;Koch et al. 2004), that is, reducing their growth efficiency with respect to the amount of light their crowns intercepted.More recently, it has been suggested that the decline in efficiency is due to respiratory costs increasing progressively with tree size (West 2020).In this theory, it is suggested that greater respiratory costs are incurred in processes such as construction of more complex anatomy of newly developing leaves, changes that ameliorate photosynthetic limitations of the canopy due to water stress.
These physiological phenomena may provide an answer to the objective of this study to explain how varying site productivities affect individual tree maximum stem growth rates (Fig. 4).Several studies, with different tree species in various parts of the world, have developed models relating tree height to stem diameter across tree populations (Curtis 1967;Krumland and Wensel 1988;Hanus et al. 1999;Calama and Montero 2004;Dorado et al. 2006;Crecente-Campo et al. 2010;Gómez-García et al. 2014;Sharma et al. 2016;Zhang et al. 2019;Raptis et al. 2021;Lin et al. 2022).They have found that trees of the same diameter are taller in stands with greater dominant heights (defined as the average height of a restricted number of the tallest or largest diameter trees in a stand).However, heights and diameters will increase progressively with age.Consider two stands of the same age but growing on sites of different productivity.The dominant height on the more productive site will be greater than that on the less productive site.This is the principle on which is based the concept of using site index (dominant height of trees at an arbitrarily chosen age, West 2015) as a measure of site productivity.Trees of the same diameter on the two sites will then be taller on the more productive site and hence, more likely to show a decline in growth rate in consequence of the physiological constraints they are suffering.Thus, it will be trees of larger diameters beyond which maximum growth rates will start to decline on the less productive site.
Teck and Hilt (1991) developed functions to relate maximum individual stem growth rates to stem size and stand site index and used growth data collected from many forest plots over periods of many years.They used data from 28 softwood  2. The line (----) shows how, in any stand, the basal area of the tree that shows the greatest growth rate of all trees in that stand tends to decline with increasing site productivity and hardwood species that occur in forests across northeastern USA.Similar work by Pretzsch and Biber (2010) considered two softwood and two hardwood species from forests of Bavaria, Germany.The function used by Teck and Hilt was: whilst Pretzsch and Biber used: where, ΔB is individual tree basal area growth rate, ΔD is diameter growth rate, D is stem diameter at the start of a growth period, S is stand site index and a, b, c and d are parameters.In both cases, the functions they fitted predicted that the overall maximum growth rate increased with increasing site index, as was case in the present study.However, Eq. 5 simply predicts maximum basal area growth rate to increase with tree size and does not predict that it will reach some overall maximum at a specific basal area.Equation 6 does allow for the determination of a size at which maximum growth rate is an overall maximum, but that size is the same for any site index.Neither of these models yielded results similar to those of Fig. 4, results that appear to show a behaviour consistent with the known physiological behaviour of tall trees.
Model (3) developed here may be used as the first stage of a more general model to predict individual basal area growth rates in relation to tree size in E. pilularis forest in its natural range of sub-tropical eastern Australia.The next stage in the development of such a model would be to determine the reductions from the maximum possible growth rates due to the effects of inter-tree competition for growth resources where trees are growing sufficiently close to each other and to determine the likelihood of mortality (Belcher et al. 1982;Wensel et al. 1987;Teck and Hilt 1991;Lessard et al. 2001;Bragg 2001;Pretzsch and Biber 2010;West 2022).Further work would be required to extend the model to predict variables of interest to forest management such as tree height, biomass or wood volume.This would often use allometric relationships between these characteristics and stem diameter or basal area (West 1979(West , 1982(West , 2015;;Paul et al. 2016).were made.The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material.If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder.To view a copy of this licence, visit http:// creat iveco mmons.org/ licen ses/ by/4.0/.
(1) ΔB = αB − B Table 1 Details of the 10 class maxima data subsets into which the data from 35,876 individual trees were splitThe plot site index ranges used to split the data into each subset are shown, together with the number of trees in each subset.The average site index is shown of the plots from which came the 50 (at most) fastest growing individuals chosen from 50 equal sized basal area classes in each site index class.The parameter estimates determined for Model (1) in each site index class are shown Site index range in class subset (m) Number of individuals in subset Average site index of selected plots (m) Parameter estimates, Model (1)

Fig. 1
Fig.1Scatter plot of stand stocking densities against stand basal areas at the start of each of the 1,042 measurement periods of the 96 plots (•), together with data from the 460 observations from plots that contained trees with maximum growth rates and that were included in the 10 class maxima data subsets (○)

Fig. 2
Fig. 2 Scatter plots (•) of individual tree basal area growth rates against tree basal areas, at the start of the growth period concerned, of the data included in four examples of the class maxima data sets constructed here from plots with site indices in the ranges a: 31-40 m, b: 26-28.5 m, c: 22-23 m and d: 13-17 m.The short, dashed lines (----) show the fit to the data in each case of the firststage Model (1), fitted with nonlinear regression, with parameter

Fig. 3
Fig.3Scatter plots (•) of the relationships between the parameter estimates of Model (1) (Table1) and the variable with which each appeared to be most closely related.The lines show the ordinary least-squares fit to the data of a: Model (2a), b: Model (2b) and c:Model (2c)

Fig. 4
Fig. 4 Predictions ( ______ ) for five stands, with site indices as indicated, of maximum basal area growth rates from Model (3) using the parameter estimates in Table2.The line (----) shows how, in any stand, the basal area of the tree that shows the greatest growth rate of all trees in that stand tends to decline with increasing site productivity

Acknowledgements
I am grateful to my colleague, Dr D.A. Ratkowsky, for valuable discussions related to methods to fit the model developed here.Journal reviewers offered a number of useful suggestions to improve the manuscript.Funding Open Access funding enabled and organized by CAUL and its Member Institutions.Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes a + bS)D c exp(−dD)

Table 2
Parameter estimates of Model (3) and their standard errors