- Email: [email protected]
Predicting the height–diameter pattern of planted Pinus kesiya stands in Zambia and Zimbabwe
Forest Ecology and Management 175 (2003) 355±366
Predicting the height±diameter pattern of planted Pinus kesiya stands in Zambia and Zimbabwe Kalle EerikaÈinen* Faculty of Forestry, University of Joensuu, P.O. Box 111, FIN-80101 Joensuu, Finland Received 11 May 2001; received in revised form 23 January 2002; accepted 15 April 2002
Abstract The aim of the study was to de®ne and analyse a prediction model for the relationship between tree height and diameter of Pinus kesiya (Royle ex Gordon) grown in the forest plantations of Zambia and Zimbabwe. The modelling was based on an equation of the power type. The development of the height±diameter pattern was assumed to be dependent on both tree size and stand characteristics. The tree diameter at breast height, stand age, dominant height and diameter were used as predictors of the tree height. Due to the spatially hierarchical (stands, trees) and temporal (measurement occasions) correlation structures in the data, the basic assumption about noncorrelated error terms did not hold. Therefore, the generalised least squares method was used in the parameter estimation of the random parameter model. In model applications, the development of the height±diameter pattern is predicted with the parameters of the ®xed model part. The variance estimates for the random stand, measurement and tree effects can be used in model calibration by applying the standard linear prediction theory. These results lead to the recommendation that the power type of model be used as the predictor of the height±diameter relationship of P. kesiya stands in the tree plantations of southeastern Africa. # 2002 Elsevier Science B.V. All rights reserved. Keywords: Power equation; Random effects; Random parameter model
1. Introduction Growth and yield simulators are generally based on statistical prediction models for the temporal development of a stand's dimensions and structure. In yield studies, stand volumes are normally calculated from tree volumes, using stem volume prediction models and trees derived from the diameter and height distributions. Horizontal stand dimensions, used in forest mensuration and growth modelling, are obtained from tree diameters. Examples of these models are: basal *
Tel.: 358-13-251-4521; fax: 358-13-251-4444. E-mail address: [email protected] (K. EerikaÈinen).
area prediction and projection models, theoretical diameter and basal area frequency distributions and, in the case of individual tree growth models, increment models for tree diameters or basal areas. Height models describe the relationship between horizontal and vertical tree and stand dimensions, i.e. tree diameters and heights and the stand characteristics derived from them. Height models can be classi®ed according to their methodological basis as follows: (1) height development models for the stand median or mean diameter trees (see e.g., Kilkki and Siitonen, 1975; Pukkala and Pohjonen, 1993); (2) static models for the relationship between tree height and diameter (e.g., Kilkki and Siitonen,
0378-1127/02/$ ± see front matter # 2002 Elsevier Science B.V. All rights reserved. PII: S 0 3 7 8 - 1 1 2 7 ( 0 2 ) 0 0 1 3 8 - X
356
K. EerikaÈinen / Forest Ecology and Management 175 (2003) 355±366
1975; Lappi, 1991, 1997; Matney and Sullivan, 1982; Flewelling and de Jong, 1994; HoÈkkaÈ, 1997); (3) height increment models of individual trees (e.g., NyyssoÈnen and MielikaÈinen, 1978; Hynynen, 1995; Hasenauer et al., 1998; Huang and Titus, 1999). In this study, the analysis was based on the height models of the second category. In the modelling of the height±diameter relationship an important factor is the temporal development of this relationship. There are generally two major elements characterising the temporal development of height±diameter curves: (1) the form development; (2) the asymptotic development of the height curve (e.g., Loetsch et al., 1973; Stage, 1975). In the case of evenaged forests, the form development usually means that the height±diameter curves are steeper in young stands and become more even towards the end of the rotation. Thus the height curve should have a steeper slope for small trees than for larger ones (Stage, 1975). Asymptotic development means that the asymptotic maximum of the curve increases as a function of age in spite of the fact that the height±diameter model may not have a speci®c asymptote parameter. It is possible to relate the asymptotic development of a height± diameter model to the development of stand dominant height. This is justi®ed because the stand dominant height is often a reliable and measurable indicator of site quality and it can therefore, be used as an indicator of the growth potential of tree heights and diameters. These factors are to be taken into account in the selection of model types and suitable parameterisations of height±diameter curves. There are several model forms available for the modelling of the height±diameter relationship. Classical model forms applied in many forestry studies are, e.g., NaÈslund's (1937) height model (e.g., Kilkki and Siitonen, 1975; Pukkala et al., 1990; Pukkala and Pohjonen, 1993) and different modi®cations of a Schumacher (1939) type of growth model (e.g., Lappi, 1991; Pukkala et al., 1998), which are basically nonlinear functions but can mathematically be transformed into a linear form. The height model used by Flewelling and de Jong (1994) is also nonlinear. The same exponential model type was also applied, after ®xing a nonlinear model parameter, by Lappi
(1997) for the longitudinal analysis of the height and diameter relationship. The height modelling approach of HoÈkkaÈ (1997) resembles and is related to Lappi's (1997) presentation. The earlier height±diameter model for Pinus kesiya (Royle ex Gordon) by Miina et al. (1999) formed the basis for this current study. However, the study of Miina et al. (1999) was based on data collected from only one tree plantation, namely Chati in Zambia. Therefore, the height±diameter model is a `local model'. The preliminary analyses of this study showed that the composition and transformations of the independent variables of Miina et al. (1999) were not adequate for a more general model. The basic nonlinear parameterisation of the earlier model was not directly usable for the modelling of a linear height± diameter model. In addition, Miina et al. (1999) could have been included more levels of random effects in their model. The aim of this study was to develop a linear random parameter model for the relationship between tree height and diameter in managed P. kesiya stands in the forest plantations of Zimbabwe and Zambia. A mixed modelling approach was emphasised due to the spatial and temporal correlation structures of the data. For the sake of a simpli®ed and controllable modelling with analyses of random effects it was required that the height±diameter model can be estimated using linear estimators. The aim was a common height prediction model that could be used: (1) in the yield simulators based on the predicted development of the stand dominant height and diameter distribution; (2) in the basic calculations of stand characteristics of the inventories in P. kesiya plantations of southeastern Africa. Thus, a special attention was paid to the parametric determination of: (a) the temporal model form; (b) the asymptotic development of the height± diameter curve, both of which were related to the biologically reasoned predictors and parametric model formulations. 2. Material The sample tree data on P. kesiya were measured from permanent growth trials in Zambia and Zimbabwe. The number and timing of successive measurement occasions and the design of sample plots
K. EerikaÈinen / Forest Ecology and Management 175 (2003) 355±366
varied between countries and plantations. The number of measured height sample trees and the criteria used in their selection varied also because of different sample plot designs. Only the measurement occasions with a minimum of 10 assessed representative trees were used in the analysis in order to achieve a reliable determination of the height pattern for every measurement occasion. The Zambian data were measured in the forest plantations in the Copperbelt region (Table 1). All Zambian permanent sample plots (PSP) (Control Plan 5/3/1, 1969) were in managed stands. The date of thinnings was known only in a few cases, and the thinning intensity was never known. Assessed tree characteristics were: diameter (mm) at breast height of all living trees and height (dm) of height sample trees, which had been selected using the rule of one tree per each 1 cm diameter class (Control Plan 5/3/1, 1969). In addition, diameters and heights of the 100 largest diameter (at breast height) trees per hectare were measured in order to determine the stand dominant diameter and height. The stand dominant diameter and height were calculated as the arithmetic means of diameters and heights of the 100 dominant trees, which is the practice in southern Africa (see also, van Laar and AkcËa, 1997, p. 247). The Zambian data was divided into the modelling and test data (Table 1). The modelling data were used for the estimation of model parameters whereas the test data were for the model validation. Naturally, the reliability of the model was also tested in the modelling data. In the Zambian data there were some stands or compartments that included more than one sample plot. Because the variation between heights and diameters within stands was low only one sample plot per stand was selected into the modelling data and the other plots were moved to the test data. Two different plantations were represented in the Zimbabwean growth data (Table 1), measured from correlated curve trend (CCT) experiments (Mushove, 1992) of the Nelder type (Nelder, 1962). Growth trials from Martin and Marondera were measured 9±13 times. Observations from these plantations were included in the modelling data (Table 1) to improve the temporal coverage of the data. Assessed tree characteristics were: diameter (mm) at breast height and height (dm) of all living trees. The health status and the technical quality of all living trees were also
357
recorded in every measurement. Only healthy trees were used in this study. When the tree and stand variables are examined the countrywise modelling data sets are rather similar (Table 1). The means of stand and tree characteristics are greater in the Zambian data than in the Zimbabwean data because of the higher number of old stands in the Zambian data (Table 1). 3. Methods 3.1. Model development In this study the development of the height±diameter model begins with the following simple power function, also known as the allometric or multiplicative model: y axb ;
(1)
where y is the dependent variable (e.g., total tree height), x the independent variable (e.g., tree diameter at breast height), and a and b are model parameters. In the case of this height±diameter model, we can assume that parameters a and b are functions of tree and stand characteristics. The functional determination of the model components a and b was based on the height increment model of Hynynen (1995) and the height±diameter model of Miina et al. (1999). These are, however, nonlinear model parameterisations. Several combinations and transformations of predictors were tested when selecting the ®xed model part of the model. The selection criteria used were the signi®cance and logical signs of the regression coef®cients. The signi®cance test of estimated model parameters tested whether the true value of the parameter is zero or not. The test value was calculated by dividing the parameter estimate by its estimated standard error and comparing this ratio to the t-distribution with n±pa degrees of freedom (n is the number of observations, and pa is the number of model parameters); if the absolute test value was greater than 1.96, the parameter was treated as signi®cant (p-value < 0:05). In modelling, one aim was to minimise the estimate of the mean square error of the model. In addition, the model residuals were analysed. In this study, the following nonlinear parameterisation
358
Data set
Plantation
Number of stands
Number of Altitude observations (m)
Rainfall (mm)
Coordinates 8S
8E
hkji (m)
dkji (cm)
Hdom kj (m)
Ddom kj (cm)
Tkj (years)
Nkj
Modelling Zambia
Chati Dola hill Lamba
94 147 10
5505 9048 558
1270 1230 1270
1300 1270 1300
12.8 13.0 12.9
27.8 28.5 27.8
Minimum Maximum Mean Median S.D.
3.0 34.2 17.8 18.0 5.371
3.0 52.8 22.0 22.0 6.309
7.2 33.1 19.2 19.4 4.965
14.3 43.7 27.7 28.1 5.051
5.6 32.1 13.5 12.8 4.763
2.0 6.0 3.9 4.0 1.054
Modelling Zimbabwe Sub total
Marondera Martin
7 7 265
3426 2758 21295
1250 1290
950 1050
18.1 19.5
31.4 32.6
Minimum Maximum Mean Median S.D.
1.4 33.0 11.4 11.6 6.464
0.8 52.3 16.7 15.9 9.251
2.4 30.4 12.6 13.4 6.500
3.0 44.3 20.8 22.9 9.735
2.2 25.8 8.7 8.7 5.110
9.0 13.0 11.2 10.5 1.528
Test Zambia
Chati Dola hill Lamba
8 25 3 36 301
383 1667 177 2227 23522
Minimum Maximum Mean Median S.D.
6.1 32.0 18.8 18.5 4.704
5.3 47.7 23.7 23.6 5.984
11.2 30.8 20.0 19.5 4.270
18.2 40.6 28.9 28.1 4.516
7.7 30.1 14.3 12.8 4.612
2.0 6.0 4.1 4.0 1.190
Sub total Total
a hkji is the total tree height, dkji the tree diameter at breast height, Hdom kj the dominant height, Ddom kj the dominant diameter, Tkj the stand age, Nkj the number of measurement occasions per sample plot, i.e. stand, and S.D. is the standard deviation of the variable.
K. EerikaÈinen / Forest Ecology and Management 175 (2003) 355±366
Table 1 Site characteristics for plantations, numbers of stands and individual height observations (hkji) for the modelling and test data of this study, and variation of tree and stand characteristics in the study materiala
K. EerikaÈinen / Forest Ecology and Management 175 (2003) 355±366
of the height±diameter model was found to provide the best ®t: dkji bkji 0 ekji ; hkji akj Ddom kj c1 in which akj C00 Hdom kj ; dkji bkji c2 c3 ln Ddom kj Hdom kj c4 ln Tkj c5 ; Tkj
(2)
where hkji is the total tree height of tree i at the measurement occasion j of stand k (m); dkji the tree diameter over bark at breast height (cm); Ddom kj the dominant diameter (cm); Hdom kj the dominant height (m); Tkj the stand age since planting (years); C00 ; c1 ; . . . ; c5 are regression coef®cients to be estimated; and e0kji is the random error term of the model. Eq. (2) can be linearised as follows: dkji ln hkji c0 c1 ln Hdom kj c2 ln Ddom kj 2 dkji dkji c4 ln c3 ln ln Tkj Ddom kj Ddom kj dkji Hdom kj c5 ln (3) e0 kji ; Ddom kj Tkj where c0 ln C00 , and e0 kji ln e0kji . Because of the spatially hierarchical (stands and trees) and temporal (measurement occasions) correlation structures of the data the basic assumption about noncorrelated residuals did not hold. Therefore, a random parameter modelling approach was applied, i.e. random effects were taken into consideration in the model formulation and estimation of ®xed and random model parameters (e.g., Lappi, 1986, 1991, 1993, 1997; Lappi and Bailey, 1988; Lappi and Malinen, 1994). The timing of stand measurements was not ®xed and it varied even within plantations. Therefore, time effects were treated as stand speci®c. On the other hand, individual tree observations of the data are structurally cross-classi®ed (Goldstein, 1995; Goldstein et al., 1998) by: (a) the spatially hierarchical stand and tree units; (b) the occasions of measurements that they pertain to. For example, a single tree may be measured several times and an observation of a tree may have been made on any of several measurement occasions.
359
In addition, there are at least 10 trees measured from every stand. Thus trees are not nested within measurements and vice versa. They are, however, nested within stands. As a result, coef®cient c4 of Eq. (3), which is related to the relative tree size and stand age, was assumed to vary randomly from stand to stand (k), from one stand speci®c measurement occasion to another (kj), and from tree to tree (ki). Therefore, the random stand, time (measurement) and tree effects were estimated for the coef®cient. The linear random parameter model for the height± diameter ratio of P. kesiya is as follows: dkji ln hkji c0 c1 ln Hdom kj c2 ln Ddom kj 2 dkji dkji c3 ln c4 kji ln ln Tkj Ddom kj Ddom kj dkji Hdom kj c5 ln e0 kji ; Ddom kj Tkj in which c4 kji c4 us4k um4kj ut4ki ;
(4)
where us4k is the random coef®cient for stand effects, um4kj the random coef®cient for time effects, ut4ki the random coef®cient for tree effects, and e0 kji is the random error term of the model. In addition, k 1; . . . ; s; j 1; . . . ; mk ; and i 1; . . . ; tk . Generally, it is assumed that cov us4k ; us4k0 cov um4kj ; um4kj0 cov ut4ki ; ut4ki0 0, when k 6 k0 ; kj 6 kj0 ; ki 6 ki0 . The estimation of random stand, time and tree effects is based on the following assumptions: E us4k E um4kj E ut4ki 0, and var us4k s2us4 , var um4kj s2um4 and var ut4ki s2ut4 . For the random error term of the model the assumptions are as follows: E e0 kji 0 and var e0 kji s2e0 . Furthermore, random effects at different correlation levels are assumed to be independent of each other, i.e. cov us4k ; um4k0 j cov us4k ; ut4k0 i cov um4kj ; ut4k0 i ; . . . ; cov ut4ki ; e0 k0 ji0 0 for all k, k0 , j, j0 , i and i0 . The parameters of the height±diameter model were estimated using the MLwiN software for multilevel analysis of linear regression models with random parameters (Goldstein et al., 1998). The method used was the iterative generalised least squares (IGLS), which provides the maximum likelihood estimates under multivariate normality (see, Goldstein, 1986, 1989).
K. EerikaÈinen / Forest Ecology and Management 175 (2003) 355±366
360
Parameters of the height±diameter model (Eq. (4)) were estimated using a transformed, i.e. linearised, model form. Thus, a bias correction is needed for the model to be added to the estimated or predicted logarithmic height before its back-transformation into the anti-logarithmic scale (see Baskerville, 1972; Flewelling and Pienaar, 1981). For the model of this study, the bias correction is half of the estimated error variance of the model, as suggested by Baskerville (1972). The estimated error variance of the model (Eq. (4)) is proportional to the tree size according to the following formula (see e.g., Lappi, 1991; HoÈkkaÈ, 1997): c c var ln hkji
s2us4
ln hkji s2um4
s2ut4
2 dkji ln ln Tkj s2e0 : Ddom kj (5)
The bias corrected height±diameter prediction model can now be expressed as follows: ^ c c hkji exp c ln hkji 12 var ln hkji
ln hkji :
(6)
3.2. Calculation of reliability ®gures In the reliability tests the means (Bias, m; Bias%, %) and standard deviations (RMSE, m; RMSE%, %) of residuals (observed estimated), and the estimates of the standard errors of residual means (S.E., m) were calculated using the following formulas: Xn hl ^ hl ; (7) Bias l1 n 100 Bias Bias% Pn ; (8) ^ l1 hl =n s Pn ^ hl 2 l1 hl ; (9) RMSE n 100 RMSE RMSE% Pn ; ^ l1 hl =n q Pn ^ hl 2 = n 1 l1 hl p S:E: n s Pn ^ hl 2 l1 hl ; n2 n
(10)
(11)
where hl is the observed value, and ^hl the estimated (modelling data) or predicted (test data) value of a tree height (m). 4. Results The ®xed parameter estimates of Eq. (4) are statistically signi®cant (Table 2). In addition, their signs are logical; the term ln dkji =Ddom kj 2 with a negative coef®cient ensures that the height±diameter relationship does not increase without bound. The variance estimates for the random stand, time related measurement and tree effects are also signi®cant (Table 2). The means (Eq. (7)) and standard deviations (Eq. (9)) of residuals, and the estimates of standard errors (Eq. (11)) of residual means for the bias corrected height±diameter model (6) were calculated for 3 cm diameter classes (Fig. 1a1 and a2), 2 m height classes of estimated tree heights (Fig. 1b1 and b2) and 2-year age classes (Fig. 1c1 and c2) in the modelling (Fig. 1a1, b1 and c1) and test (Fig. 1a2, b2 and c2) data. When the reliability ®gures are inspected in the modelling data by diameter (Fig. 1a1) and height (Fig. 1b1) classes, it can be seen that the height±diameter prediction model gives almost unbiased estimates for trees of all sizes. However, the model tends to overestimate the heights of the largest trees. The height predictions are slightly better for the test data than for the modelling data, i.e. the means of the residuals are
Table 2 Statistics for the height±diameter model (Eq. (4))a Parameter c0 c1 c2 c3 c4 c5 s2us4 s2um4 s2ut4 s2e0
Estimate 0.038431 1.014449 0.228008 0.029235 0.120584 0.294715 0.000726 0.001853 0.005930 0.005109
Standard error 0.004691 0.001664 0.040968 0.008730 0.007219 0.023196 0.000149 0.000121 0.000187 0.000059
a 2 sus4 is the variance estimate for random stand effects, s2um4 the variance estimate for random measurement effects, s2ut4 the variance estimate for random tree effects, and s2e0 the variance estimate for random errors of the model.
K. EerikaÈinen / Forest Ecology and Management 175 (2003) 355±366
361
Fig. 1. Residual means (-^-, Eq. (7)) and deviations ( , Eq. (9)), and estimates of standard errors of residual means (- - -, Eq. (11)) in the modelling (a1, b1, c1) and test (a2, b2, c2) data with a respect to tree diameter at breast height in 3 cm diameter classes (a1 and a2), estimated height in 2 m height classes (b1 and b2), and stand age in 2-year age classes (c1 and c2). Residuals are obtained for the anti-logarithmic, i.e. back-transformed and bias corrected form of the height±diameter model (Eq. (6)).
smaller, with the exception of the smallest trees (Fig. 1a2 and b2). The means and deviations of residuals analysed against the stand age in the modelling data reveal that the height±diameter model gives unbiased height estimates for trees up to 22 years in age but overestimates the heights beyond 22 years (Fig. 1c1). The height±diameter model also gives slight overpredictions for the trees of the oldest stands in the test data (Fig. 1c2).
Further, reliability ®gures (Eqs. (7)±(10)) were calculated for all observations in the Zambian modelling and test data and also for the Zimbabwean modelling data. Absolute and relative biases for the Zambian modelling and test data were zero. The RMSEs were 1.53 m (8.6%) for the modelling data and 1.42 m (7.5%) for the test data. For the Zimbabwean modelling data the bias was 0.1 m ( 0.4%) and the RMSE was 1.19 m (10.4%).
362
K. EerikaÈinen / Forest Ecology and Management 175 (2003) 355±366
Fig. 2. Relationship between the predicted height±diameter curve (Ð), which is obtained using Eq. (6), and measured dominant height± diameter points (&) (a1 and b1), and the predicted height±diameter curve (Ð) and observed height±diameter points (a2 and b2). a1 and a2 are obtained for the PSP from the tree plantation of Lamba in Zambia; the observed tree diameters and heights are measured at the age of 7.7 (^), 13.8 ( ) and 30.1 () years, respectively. b1 and b2 are obtained for the CCT-treatment from the tree plantation of Martin in Zimbabwe; the observed tree diameters and heights are measured at the ages of 2.2 (^), 4.2 ( ), 8.5 (), 12.4 (^), 18.4 (&) and 25.8 (~) years, respectively.
The development of the predicted height±diameter pattern was illustrated with a PSP-trial from Zambia (Fig. 2a1 and a2) and a CCT-treatment from Zimbabwe (Fig. 2b1 and b2). Fig. 2 shows that both the temporal form and asymptotic development of the height±diameter curve is attained with the model presented. The predicted height±diameter curves are steeper for the earliest stages of stands when compared to the more
even curves obtained for the later measurement occasions. 5. Discussion The aim of this study was to develop a prediction model for the height±diameter pattern of P. kesiya to
K. EerikaÈinen / Forest Ecology and Management 175 (2003) 355±366
be applied in yield simulators and basic calculations of forest inventories. Thus the model predictors used are variables typically measured in forest inventories and also commonly used in yield simulators. Special attention was paid to the modelling of the temporal development of the height±diameter pattern, for which the re-measured sample plot data were very appropriate. According to the results, the presented model predicts both the temporal form and asymptotic development of the height±diameter relationship (see Fig. 2). Steeper height±diameter curves, obtained for younger stands, develop towards more even curves of older stands. In addition, the predicted curves follow the development of stand dominant trees (see Fig. 2). These are the properties required for height± diameter models, which are used in yield simulators based on the predicted development of the stand dominant height and diameter distribution. When the means (Biases) and standard deviations (RMSEs) of model residuals were inspected in diameter (Fig. 1a1 and a2), height (Fig. 1b1 and b2) and age (Fig. 1c1 and c2) classes of the modelling (Fig. 1a1, b1 and c1) and test (Fig. 1a2, b2 and c2) data, it was observed that, even if the variations of classi®ed residuals (RMSEs) are relatively high, the model gives almost unbiased estimates for trees of all sizes. However, the height±diameter model tends to overestimate the heights of the largest trees and gives slight overpredictions for trees of the oldest stands. The overpredictions in the modelling and test data are most probably caused by an insuf®cient number of observations from old stands. In most studies on growth modelling the temporal tree and stand development is linked with the development of dominant height (e.g., Borders, 1989; SaramaÈki, 1992; Hynynen, 1995). This is justi®ed because the dominant height is a measurable stand characteristic and it indicates the site quality in terms of the stand growth rate and yield capacity. It can therefore be used as a predictor of other tree and stand variables needed for the stand growth and yield predictions (e.g., Clutter et al., 1983). With most tree species there is a strong allometric relationship between tree diameter, height and stem volume (e.g., Crow and Schlaegel, 1988). Therefore, prediction errors of stem volumes can be decreased signi®cantly by using the diameter±height ratio as the determinant of the volume (e.g., Laasasenaho, 1982).
363
In fact, the prediction of tree heights is an essential stage when tree taper models are used to obtain stand level volume yields of different timber assortments. However, it is important to consider the risk of using an explanatory variable that contain prediction error, which is obvious in the applications of standard stem volume equations when predicted tree heights are used as explanatory variables of volume models (see e.g., Gregoire and Williams, 1992). There are some earlier studies in which the modelling of the height±diameter relationship has been based on random parameter models, i.e. mixed models, with random and ®xed parameters (Lappi, 1997; HoÈkkaÈ, 1997; Miina et al., 1999). In the estimation of tree height curves, the development of the height± diameter relationship and its linking with the growth of dominant trees is easier to obtain using conditioned forms of nonlinear models (Hynynen, 1995; Miina et al., 1999). Thus, if linear regression estimators and model types are not required, it is possible to use more complex nonlinear model forms and conditions for parameters in the modelling. Parameters are then estimated using nonlinear regression estimators. However, the main aim of this study was to analyse and de®ne a height±diameter model that can be estimated using linear estimators and used in growth and yield simulators based on the development of dominant trees. If the data do not contain spatial or temporal correlation structures, i.e. observations are independent and noncorrelated, the ordinary least squares (OLS) can be used in the estimation of model parameters. However, if these assumptions are false then it is more recommendable to use GLS as an estimator (e.g., Lappi, 1986, 1991, 1997). Nevertheless, the selection of the regression estimator depends on the model forms chosen and, especially, assumptions about error terms and model structures. According to the signi®cance of the variance estimates for the random stand and tree effects of the height±diameter model, it was very evident that observations were correlated spatially. There were also temporal correlation structures caused by re-measurements. Because the basic assumption about noncorrelated residuals did not hold, the random parameter modelling approach with the linear type of model was applied (e.g., Lappi, 1991, 1997; Searle, 1987; Goldstein, 1995; Goldstein et al., 1998). Because of the repeated measurements of
364
K. EerikaÈinen / Forest Ecology and Management 175 (2003) 355±366
the sample plots the hierarchical data structure was broken down, resulting in the cross-classi®ed structure of the data. However, the data enabled the estimation of variances for the random stand, time and tree effects. The height±diameter model of this study is dependent on the dominant height±diameter relationship and its development. The dominant height±diameter relationship is used to scale the curves that are obtained with the model but these curves do not necessarily pass exactly through the observed dominant height±diameter point even if they approximate it, which is logical. Further, the ®xed part of the random parameter model (Eq. (4)) can be used when measured or predicted independent model variables are available. If the model presented is used, e.g. in forest inventories, and if assessed heights of sample trees and stand characteristics, i.e. independent model variables, are available, the model can be calibrated by applying the standard linear prediction theory. Then parameters for random effects are obtained using the best linear unbiased predictor (BLUP) (e.g., Searle, 1971, 1987; Henderson, 1975; Lappi, 1986, 1991; Lappi and Bailey, 1988). Special attention was paid to the selection of the independent variables of the height±diameter prediction model. In this study, the stand age, and dominant height and diameter were used to guarantee the logical temporal development of the height±diameter pattern. It is obvious that the parameters of Eq. (4) are highly correlated. Thus even small changes in parameterisations, i.e. transformations and combinations of predictors, were re¯ected directly in the signs and signi®cance of other parameters estimated. Finally, a simple model form was sought aiming at logical signs and signi®cant estimates of the eventual model parameters. In many studies the logical behaviour of the model with small trees is guaranteed by adding a constant 1.3 to the right hand side of the equation (e.g., NaÈslund, 1937; HoÈkkaÈ, 1997; Forss et al., 1998). By doing this, the model attains height values that are equal to or greater than 1.3 m in every circumstance. In the case of height±diameter equations, which can be transformed into a linear form by taking a logarithm, the constant value, i.e. 1.3, is transferred to the left hand side of equation before taking the logarithm. According
to the observations of Lappi (1991, 1997), however, ln h 1:3 is a very unstable dependent variable for small trees, which degreases estimation properties of a height±diameter equation when data contains measurements of small trees. This was also evident in this study. When the model presented (Eq. (6)) was analysed using the data of this study, it was observed to perform fairly enough with small trees even without the above standard conditioning of the model. The lowest height prediction 1.0 m was obtained for the tree with a diameter at breast height of 0.8 cm. In general, predicted heights attain values close to 1.3 m when the diameter at breast height approaches 1 cm. Finally, the moment when most of trees attain the height of 1.3 m is at the early stage of the rotation, which is a less important period to be considered in stand growth and yield predictions. In the reliability analysis of this study it was observed that the power type of height±diameter model (Eq. (6)) is a reliable and stable predictor. The model has a basic form that smoothes the temporal development of the height patterns of P. kesiya. In addition, the site effect is embedded in the model via the dominant height±diameter dependency of the presented random parameter height±diameter model. These results lead to the recommendation that the model of this study be used as a predictor of the height±diameter relationship of managed stands of P. kesiya plantations in Zambia and Zimbabwe. The height±diameter model presented is an applicable tool with forest inventories and, in addition, with the growth and yield simulators based, e.g., on the predicted development of diameter distribution and dominant height. Acknowledgements I wish to thank Professor Timo Pukkala, Dr. Juha Lappi and Dr. Jari Miina for their valuable comments on this study. The ®eld assessments and data collection in Africa was ®nanced by the Finnish Society of Forest Science and the University of Joensuu, Faculty of Forestry. Financial support by the Finnish Foresters' Foundation is also gratefully acknowledged. The following persons and institutions have supported me in the data collection and ®eld assessments: Mr. Progress M. Sekeli (Division of Forest Research, Zambia),
K. EerikaÈinen / Forest Ecology and Management 175 (2003) 355±366
Mr. Danaza Mabvurira (the Forestry Commission of Zimbabwe), Dr. Ladislaus Nshubemuki (the Tanzania Forestry Research Institute) and Dr. Jussi SaramaÈki. This study is a part of the INCO DC cooperation project Tree Seedling Production and Management of Plantation Forests. Finally, I wish to thank Dr. Lisa Lena Opas-HaÈnninen for revising my English. References Baskerville, G.L., 1972. Use of logarithmic regression in the estimation of plant biomass. Can. J. For. Res. 2, 49±53. Borders, B.E., 1989. System of equations in forest stand modelling. For. Sci. 35, 548±556. Clutter, J.L., Fortson, J.C., Pienaar, L.V., Brister, G.H., Bailey, R.L., 1983. Timber Management: A Quantitative Approach. Wiley, New York, 333 pp. Control Plan 5/3/1, 1969. Systems of measurements in commercial compartments. A System of Continuous Inventory for Pinus Species Established by Industrial Plantations. Mimeograph, Zambia Forest Department, Division of Forest Research, Zambia, 5 pp. Crow, T.R., Schlaegel, B.E., 1988. A guide to using regression equations for estimating tree biomass. North. J. Appl. For. 5, 15±22. Flewelling, J.W., de Jong, R., 1994. Considerations in simultaneous curve ®tting for repeated height±diameter measurements. Can. J. For. Res. 24, 1408±1414. Flewelling, J.W., Pienaar, L.V., 1981. Multiplicative regression with lognormal errors. For. Sci. 27, 281±289. Forss, E., Maltamo, M., SaramaÈki, J., 1998. Static stand and tree characteristics models for Acacia mangium plantations. J. Trop. For. Sci. 10, 318±336. Goldstein, H., 1986. Multilevel mixed linear model analysis using iterative generalized least squares. Biometrika 73, 43±56. Goldstein, H., 1989. Restricted unbiased iterative generalized leastsquares estimation. Biometrika 76, 622±623. Goldstein, H., 1995. Multilevel Statistical Models, 2nd Edition. Institute of Education, University of London, 178 pp. Goldstein, H., Rasbash, J., Plewis, I., Draper, D., Browne, W., Yang, M., Woodhouse, G., Healy, M., 1998. A User's Guide to MLwiN, Multilevel Models Project. Institute of Education, University of London, 140 pp. Gregoire, T.G., Williams, M., 1992. Identifying and evaluating the components of non-measurement error in the application of standard volume equations. Statistician 41, 509±518. Hasenauer, H., Monserud, R.A., Gregoire, T.G., 1998. Using simultaneous regression techniques with individual-tree growth models. For. Sci. 44, 87±95. Henderson, C.R., 1975. Best linear unbiased estimation and prediction under a selection model. Biometrics 31, 423±447. HoÈkkaÈ, H., 1997. Height±diameter curves with random intercepts and slopes for trees growing on drained peatlands. For. Ecol. Manage. 97, 63±72.
365
Huang, S., Titus, S.J., 1999. An individual tree height increment model for mixed white spruce-aspen stands in Alberta. Can. For. Ecol. Manage. 123, 41±53. Hynynen, J., 1995. Predicting the growth response to thinning for Scots pine stands using individual-tree growth models. Silva Fenn. 29, 225±246. Kilkki, P., Siitonen, M., 1975. Simulation of arti®cial stands and derivation of growing stock models from this material. Acta For. Fenn. 145, 33, (in Finnish with English summary). Laasasenaho, J., 1982. Taper curve and volume functions for pine, spruce and birch. Commun. Inst. For. Fenn. 108, 74. Lappi, J., 1986. Mixed linear models for analyzing and predicting stem form variation of Scots pine. Commun. Inst. For. Fenn. 134, 69. Lappi, J., 1991. Calibration of height and volume equations with random parameters. For. Sci. 37, 781±801. Lappi, J., 1993. Biometrical Methods for Forest Sciences. Faculty of Forestry, University of Joensuu. University Press, Joensuu, Finland. Silva Carelica 24, 182 pp. (in Finnish). Lappi, J., 1997. A longitudinal analysis of height/diameter curves. For. Sci. 43, 555±570. Lappi, J., Bailey, R.L., 1988. A height prediction model with random stand and tree parameters: an alternative to traditional site index methods. For. Sci. 34, 907±927. Lappi, J., Malinen, J., 1994. Random-parameter height/age models when stand parameters and stand age are correlated. For. Sci. 40, 715±731. Loetsch, F., ZoÈhrer, F., Haller, K.E., 1973. Forest Inventory: 2. BLV Verlagsgesellschaft, MuÈnchen, Germany, 469 pp. Matney, T.G., Sullivan, A.D., 1982. Compatible stand and stock tables for thinned and unthinned loblolly pine stands. For. Sci. 28, 161±171. Miina, J., Maltamo, M., EerikaÈinen, K., 1999. Modelling the growth and yield of Pinus kesiya in Chati, Zambia. In: Pukkala, T., EerikaÈinen, K. (Eds.), Growth and Yield Modelling of Tree Plantations in South and East Africa. Research Notes 97. Faculty of Forestry, University of Joensuu, Finland, pp. 167±181. Mushove, P.T., 1992. Thirty years of correlated curve trend management trials. In: Piearce, G.D., Shaw, P. (Eds.), Proceedings of the Anniversary Seminar on Forestry Research in Zimbabwe, Mutare, The Forestry Commission, Harare, Zimbabwe, August 27±31, 1990, pp. 75±83. NaÈslund, M., 1937. SkogsfoÈrsoÈksanstaltens gallringsfoÈrsoÈk i tallskog. Rep. Swed. Inst. Exp. For. 29, 1±169, (in Swedish with German summary). Nelder, J.A., 1962. New kinds of systematic designs for spacing experiments. Biometrics 18, 283±307. NyyssoÈnen, A., MielikaÈinen, K., 1978. Estimation of stand increment. Acta For. Fenn. 163, 40, (in Finnish with English summary). Pukkala, T., Pohjonen, V., 1993. Yield of Cupressus lusitanica in Ethiopia. Silva Fenn. 27, 195±207. Pukkala, T., SaramaÈki, J., Mubita, O., 1990. Management planning system for tree plantations. A case study for Pinus kesiya in Zambia. Silva Fenn. 24, 171±180. Pukkala, T., Miina, J., Kurttila, M., KolstroÈm, T., 1998. A spatial yield model for optimizing the thinning regime of mixed stands
366
K. EerikaÈinen / Forest Ecology and Management 175 (2003) 355±366
of Pinus sylvestris and Picea abies. Scand. J. For. Res. 13, 31± 42. SaramaÈki, J., 1992. A growth and yield prediction model of Pinus kesiya (Royle ex Gordon) in Zambia. Acta For. Fenn. 230, 68. Schumacher, F.X., 1939. A new growth curve and its application to timber-yield studies. J. For. 37, 819±820.
Searle, S.R., 1971. Linear Models. Wiley, New York, 532 pp. Searle, S.R., 1987. Linear Models for Unbalanced Data. Wiley, New York, 536 pp. Stage, A.R., 1975. Prediction of height increment for models of forest growth. USDA For. Serv. Res. Paper INT-164, 20 pp. van Laar, A., AkcËa, A., 1997. Forest Mensuration. Cuvillier Verlag, GoÈttingen, 418 pp.
Predicting the height±diameter pattern of planted Pinus kesiya stands in Zambia and Zimbabwe Kalle EerikaÈinen* Faculty of Forestry, University of Joensuu, P.O. Box 111, FIN-80101 Joensuu, Finland Received 11 May 2001; received in revised form 23 January 2002; accepted 15 April 2002
Abstract The aim of the study was to de®ne and analyse a prediction model for the relationship between tree height and diameter of Pinus kesiya (Royle ex Gordon) grown in the forest plantations of Zambia and Zimbabwe. The modelling was based on an equation of the power type. The development of the height±diameter pattern was assumed to be dependent on both tree size and stand characteristics. The tree diameter at breast height, stand age, dominant height and diameter were used as predictors of the tree height. Due to the spatially hierarchical (stands, trees) and temporal (measurement occasions) correlation structures in the data, the basic assumption about noncorrelated error terms did not hold. Therefore, the generalised least squares method was used in the parameter estimation of the random parameter model. In model applications, the development of the height±diameter pattern is predicted with the parameters of the ®xed model part. The variance estimates for the random stand, measurement and tree effects can be used in model calibration by applying the standard linear prediction theory. These results lead to the recommendation that the power type of model be used as the predictor of the height±diameter relationship of P. kesiya stands in the tree plantations of southeastern Africa. # 2002 Elsevier Science B.V. All rights reserved. Keywords: Power equation; Random effects; Random parameter model
1. Introduction Growth and yield simulators are generally based on statistical prediction models for the temporal development of a stand's dimensions and structure. In yield studies, stand volumes are normally calculated from tree volumes, using stem volume prediction models and trees derived from the diameter and height distributions. Horizontal stand dimensions, used in forest mensuration and growth modelling, are obtained from tree diameters. Examples of these models are: basal *
Tel.: 358-13-251-4521; fax: 358-13-251-4444. E-mail address: [email protected] (K. EerikaÈinen).
area prediction and projection models, theoretical diameter and basal area frequency distributions and, in the case of individual tree growth models, increment models for tree diameters or basal areas. Height models describe the relationship between horizontal and vertical tree and stand dimensions, i.e. tree diameters and heights and the stand characteristics derived from them. Height models can be classi®ed according to their methodological basis as follows: (1) height development models for the stand median or mean diameter trees (see e.g., Kilkki and Siitonen, 1975; Pukkala and Pohjonen, 1993); (2) static models for the relationship between tree height and diameter (e.g., Kilkki and Siitonen,
0378-1127/02/$ ± see front matter # 2002 Elsevier Science B.V. All rights reserved. PII: S 0 3 7 8 - 1 1 2 7 ( 0 2 ) 0 0 1 3 8 - X
356
K. EerikaÈinen / Forest Ecology and Management 175 (2003) 355±366
1975; Lappi, 1991, 1997; Matney and Sullivan, 1982; Flewelling and de Jong, 1994; HoÈkkaÈ, 1997); (3) height increment models of individual trees (e.g., NyyssoÈnen and MielikaÈinen, 1978; Hynynen, 1995; Hasenauer et al., 1998; Huang and Titus, 1999). In this study, the analysis was based on the height models of the second category. In the modelling of the height±diameter relationship an important factor is the temporal development of this relationship. There are generally two major elements characterising the temporal development of height±diameter curves: (1) the form development; (2) the asymptotic development of the height curve (e.g., Loetsch et al., 1973; Stage, 1975). In the case of evenaged forests, the form development usually means that the height±diameter curves are steeper in young stands and become more even towards the end of the rotation. Thus the height curve should have a steeper slope for small trees than for larger ones (Stage, 1975). Asymptotic development means that the asymptotic maximum of the curve increases as a function of age in spite of the fact that the height±diameter model may not have a speci®c asymptote parameter. It is possible to relate the asymptotic development of a height± diameter model to the development of stand dominant height. This is justi®ed because the stand dominant height is often a reliable and measurable indicator of site quality and it can therefore, be used as an indicator of the growth potential of tree heights and diameters. These factors are to be taken into account in the selection of model types and suitable parameterisations of height±diameter curves. There are several model forms available for the modelling of the height±diameter relationship. Classical model forms applied in many forestry studies are, e.g., NaÈslund's (1937) height model (e.g., Kilkki and Siitonen, 1975; Pukkala et al., 1990; Pukkala and Pohjonen, 1993) and different modi®cations of a Schumacher (1939) type of growth model (e.g., Lappi, 1991; Pukkala et al., 1998), which are basically nonlinear functions but can mathematically be transformed into a linear form. The height model used by Flewelling and de Jong (1994) is also nonlinear. The same exponential model type was also applied, after ®xing a nonlinear model parameter, by Lappi
(1997) for the longitudinal analysis of the height and diameter relationship. The height modelling approach of HoÈkkaÈ (1997) resembles and is related to Lappi's (1997) presentation. The earlier height±diameter model for Pinus kesiya (Royle ex Gordon) by Miina et al. (1999) formed the basis for this current study. However, the study of Miina et al. (1999) was based on data collected from only one tree plantation, namely Chati in Zambia. Therefore, the height±diameter model is a `local model'. The preliminary analyses of this study showed that the composition and transformations of the independent variables of Miina et al. (1999) were not adequate for a more general model. The basic nonlinear parameterisation of the earlier model was not directly usable for the modelling of a linear height± diameter model. In addition, Miina et al. (1999) could have been included more levels of random effects in their model. The aim of this study was to develop a linear random parameter model for the relationship between tree height and diameter in managed P. kesiya stands in the forest plantations of Zimbabwe and Zambia. A mixed modelling approach was emphasised due to the spatial and temporal correlation structures of the data. For the sake of a simpli®ed and controllable modelling with analyses of random effects it was required that the height±diameter model can be estimated using linear estimators. The aim was a common height prediction model that could be used: (1) in the yield simulators based on the predicted development of the stand dominant height and diameter distribution; (2) in the basic calculations of stand characteristics of the inventories in P. kesiya plantations of southeastern Africa. Thus, a special attention was paid to the parametric determination of: (a) the temporal model form; (b) the asymptotic development of the height± diameter curve, both of which were related to the biologically reasoned predictors and parametric model formulations. 2. Material The sample tree data on P. kesiya were measured from permanent growth trials in Zambia and Zimbabwe. The number and timing of successive measurement occasions and the design of sample plots
K. EerikaÈinen / Forest Ecology and Management 175 (2003) 355±366
varied between countries and plantations. The number of measured height sample trees and the criteria used in their selection varied also because of different sample plot designs. Only the measurement occasions with a minimum of 10 assessed representative trees were used in the analysis in order to achieve a reliable determination of the height pattern for every measurement occasion. The Zambian data were measured in the forest plantations in the Copperbelt region (Table 1). All Zambian permanent sample plots (PSP) (Control Plan 5/3/1, 1969) were in managed stands. The date of thinnings was known only in a few cases, and the thinning intensity was never known. Assessed tree characteristics were: diameter (mm) at breast height of all living trees and height (dm) of height sample trees, which had been selected using the rule of one tree per each 1 cm diameter class (Control Plan 5/3/1, 1969). In addition, diameters and heights of the 100 largest diameter (at breast height) trees per hectare were measured in order to determine the stand dominant diameter and height. The stand dominant diameter and height were calculated as the arithmetic means of diameters and heights of the 100 dominant trees, which is the practice in southern Africa (see also, van Laar and AkcËa, 1997, p. 247). The Zambian data was divided into the modelling and test data (Table 1). The modelling data were used for the estimation of model parameters whereas the test data were for the model validation. Naturally, the reliability of the model was also tested in the modelling data. In the Zambian data there were some stands or compartments that included more than one sample plot. Because the variation between heights and diameters within stands was low only one sample plot per stand was selected into the modelling data and the other plots were moved to the test data. Two different plantations were represented in the Zimbabwean growth data (Table 1), measured from correlated curve trend (CCT) experiments (Mushove, 1992) of the Nelder type (Nelder, 1962). Growth trials from Martin and Marondera were measured 9±13 times. Observations from these plantations were included in the modelling data (Table 1) to improve the temporal coverage of the data. Assessed tree characteristics were: diameter (mm) at breast height and height (dm) of all living trees. The health status and the technical quality of all living trees were also
357
recorded in every measurement. Only healthy trees were used in this study. When the tree and stand variables are examined the countrywise modelling data sets are rather similar (Table 1). The means of stand and tree characteristics are greater in the Zambian data than in the Zimbabwean data because of the higher number of old stands in the Zambian data (Table 1). 3. Methods 3.1. Model development In this study the development of the height±diameter model begins with the following simple power function, also known as the allometric or multiplicative model: y axb ;
(1)
where y is the dependent variable (e.g., total tree height), x the independent variable (e.g., tree diameter at breast height), and a and b are model parameters. In the case of this height±diameter model, we can assume that parameters a and b are functions of tree and stand characteristics. The functional determination of the model components a and b was based on the height increment model of Hynynen (1995) and the height±diameter model of Miina et al. (1999). These are, however, nonlinear model parameterisations. Several combinations and transformations of predictors were tested when selecting the ®xed model part of the model. The selection criteria used were the signi®cance and logical signs of the regression coef®cients. The signi®cance test of estimated model parameters tested whether the true value of the parameter is zero or not. The test value was calculated by dividing the parameter estimate by its estimated standard error and comparing this ratio to the t-distribution with n±pa degrees of freedom (n is the number of observations, and pa is the number of model parameters); if the absolute test value was greater than 1.96, the parameter was treated as signi®cant (p-value < 0:05). In modelling, one aim was to minimise the estimate of the mean square error of the model. In addition, the model residuals were analysed. In this study, the following nonlinear parameterisation
358
Data set
Plantation
Number of stands
Number of Altitude observations (m)
Rainfall (mm)
Coordinates 8S
8E
hkji (m)
dkji (cm)
Hdom kj (m)
Ddom kj (cm)
Tkj (years)
Nkj
Modelling Zambia
Chati Dola hill Lamba
94 147 10
5505 9048 558
1270 1230 1270
1300 1270 1300
12.8 13.0 12.9
27.8 28.5 27.8
Minimum Maximum Mean Median S.D.
3.0 34.2 17.8 18.0 5.371
3.0 52.8 22.0 22.0 6.309
7.2 33.1 19.2 19.4 4.965
14.3 43.7 27.7 28.1 5.051
5.6 32.1 13.5 12.8 4.763
2.0 6.0 3.9 4.0 1.054
Modelling Zimbabwe Sub total
Marondera Martin
7 7 265
3426 2758 21295
1250 1290
950 1050
18.1 19.5
31.4 32.6
Minimum Maximum Mean Median S.D.
1.4 33.0 11.4 11.6 6.464
0.8 52.3 16.7 15.9 9.251
2.4 30.4 12.6 13.4 6.500
3.0 44.3 20.8 22.9 9.735
2.2 25.8 8.7 8.7 5.110
9.0 13.0 11.2 10.5 1.528
Test Zambia
Chati Dola hill Lamba
8 25 3 36 301
383 1667 177 2227 23522
Minimum Maximum Mean Median S.D.
6.1 32.0 18.8 18.5 4.704
5.3 47.7 23.7 23.6 5.984
11.2 30.8 20.0 19.5 4.270
18.2 40.6 28.9 28.1 4.516
7.7 30.1 14.3 12.8 4.612
2.0 6.0 4.1 4.0 1.190
Sub total Total
a hkji is the total tree height, dkji the tree diameter at breast height, Hdom kj the dominant height, Ddom kj the dominant diameter, Tkj the stand age, Nkj the number of measurement occasions per sample plot, i.e. stand, and S.D. is the standard deviation of the variable.
K. EerikaÈinen / Forest Ecology and Management 175 (2003) 355±366
Table 1 Site characteristics for plantations, numbers of stands and individual height observations (hkji) for the modelling and test data of this study, and variation of tree and stand characteristics in the study materiala
K. EerikaÈinen / Forest Ecology and Management 175 (2003) 355±366
of the height±diameter model was found to provide the best ®t: dkji bkji 0 ekji ; hkji akj Ddom kj c1 in which akj C00 Hdom kj ; dkji bkji c2 c3 ln Ddom kj Hdom kj c4 ln Tkj c5 ; Tkj
(2)
where hkji is the total tree height of tree i at the measurement occasion j of stand k (m); dkji the tree diameter over bark at breast height (cm); Ddom kj the dominant diameter (cm); Hdom kj the dominant height (m); Tkj the stand age since planting (years); C00 ; c1 ; . . . ; c5 are regression coef®cients to be estimated; and e0kji is the random error term of the model. Eq. (2) can be linearised as follows: dkji ln hkji c0 c1 ln Hdom kj c2 ln Ddom kj 2 dkji dkji c4 ln c3 ln ln Tkj Ddom kj Ddom kj dkji Hdom kj c5 ln (3) e0 kji ; Ddom kj Tkj where c0 ln C00 , and e0 kji ln e0kji . Because of the spatially hierarchical (stands and trees) and temporal (measurement occasions) correlation structures of the data the basic assumption about noncorrelated residuals did not hold. Therefore, a random parameter modelling approach was applied, i.e. random effects were taken into consideration in the model formulation and estimation of ®xed and random model parameters (e.g., Lappi, 1986, 1991, 1993, 1997; Lappi and Bailey, 1988; Lappi and Malinen, 1994). The timing of stand measurements was not ®xed and it varied even within plantations. Therefore, time effects were treated as stand speci®c. On the other hand, individual tree observations of the data are structurally cross-classi®ed (Goldstein, 1995; Goldstein et al., 1998) by: (a) the spatially hierarchical stand and tree units; (b) the occasions of measurements that they pertain to. For example, a single tree may be measured several times and an observation of a tree may have been made on any of several measurement occasions.
359
In addition, there are at least 10 trees measured from every stand. Thus trees are not nested within measurements and vice versa. They are, however, nested within stands. As a result, coef®cient c4 of Eq. (3), which is related to the relative tree size and stand age, was assumed to vary randomly from stand to stand (k), from one stand speci®c measurement occasion to another (kj), and from tree to tree (ki). Therefore, the random stand, time (measurement) and tree effects were estimated for the coef®cient. The linear random parameter model for the height± diameter ratio of P. kesiya is as follows: dkji ln hkji c0 c1 ln Hdom kj c2 ln Ddom kj 2 dkji dkji c3 ln c4 kji ln ln Tkj Ddom kj Ddom kj dkji Hdom kj c5 ln e0 kji ; Ddom kj Tkj in which c4 kji c4 us4k um4kj ut4ki ;
(4)
where us4k is the random coef®cient for stand effects, um4kj the random coef®cient for time effects, ut4ki the random coef®cient for tree effects, and e0 kji is the random error term of the model. In addition, k 1; . . . ; s; j 1; . . . ; mk ; and i 1; . . . ; tk . Generally, it is assumed that cov us4k ; us4k0 cov um4kj ; um4kj0 cov ut4ki ; ut4ki0 0, when k 6 k0 ; kj 6 kj0 ; ki 6 ki0 . The estimation of random stand, time and tree effects is based on the following assumptions: E us4k E um4kj E ut4ki 0, and var us4k s2us4 , var um4kj s2um4 and var ut4ki s2ut4 . For the random error term of the model the assumptions are as follows: E e0 kji 0 and var e0 kji s2e0 . Furthermore, random effects at different correlation levels are assumed to be independent of each other, i.e. cov us4k ; um4k0 j cov us4k ; ut4k0 i cov um4kj ; ut4k0 i ; . . . ; cov ut4ki ; e0 k0 ji0 0 for all k, k0 , j, j0 , i and i0 . The parameters of the height±diameter model were estimated using the MLwiN software for multilevel analysis of linear regression models with random parameters (Goldstein et al., 1998). The method used was the iterative generalised least squares (IGLS), which provides the maximum likelihood estimates under multivariate normality (see, Goldstein, 1986, 1989).
K. EerikaÈinen / Forest Ecology and Management 175 (2003) 355±366
360
Parameters of the height±diameter model (Eq. (4)) were estimated using a transformed, i.e. linearised, model form. Thus, a bias correction is needed for the model to be added to the estimated or predicted logarithmic height before its back-transformation into the anti-logarithmic scale (see Baskerville, 1972; Flewelling and Pienaar, 1981). For the model of this study, the bias correction is half of the estimated error variance of the model, as suggested by Baskerville (1972). The estimated error variance of the model (Eq. (4)) is proportional to the tree size according to the following formula (see e.g., Lappi, 1991; HoÈkkaÈ, 1997): c c var ln hkji
s2us4
ln hkji s2um4
s2ut4
2 dkji ln ln Tkj s2e0 : Ddom kj (5)
The bias corrected height±diameter prediction model can now be expressed as follows: ^ c c hkji exp c ln hkji 12 var ln hkji
ln hkji :
(6)
3.2. Calculation of reliability ®gures In the reliability tests the means (Bias, m; Bias%, %) and standard deviations (RMSE, m; RMSE%, %) of residuals (observed estimated), and the estimates of the standard errors of residual means (S.E., m) were calculated using the following formulas: Xn hl ^ hl ; (7) Bias l1 n 100 Bias Bias% Pn ; (8) ^ l1 hl =n s Pn ^ hl 2 l1 hl ; (9) RMSE n 100 RMSE RMSE% Pn ; ^ l1 hl =n q Pn ^ hl 2 = n 1 l1 hl p S:E: n s Pn ^ hl 2 l1 hl ; n2 n
(10)
(11)
where hl is the observed value, and ^hl the estimated (modelling data) or predicted (test data) value of a tree height (m). 4. Results The ®xed parameter estimates of Eq. (4) are statistically signi®cant (Table 2). In addition, their signs are logical; the term ln dkji =Ddom kj 2 with a negative coef®cient ensures that the height±diameter relationship does not increase without bound. The variance estimates for the random stand, time related measurement and tree effects are also signi®cant (Table 2). The means (Eq. (7)) and standard deviations (Eq. (9)) of residuals, and the estimates of standard errors (Eq. (11)) of residual means for the bias corrected height±diameter model (6) were calculated for 3 cm diameter classes (Fig. 1a1 and a2), 2 m height classes of estimated tree heights (Fig. 1b1 and b2) and 2-year age classes (Fig. 1c1 and c2) in the modelling (Fig. 1a1, b1 and c1) and test (Fig. 1a2, b2 and c2) data. When the reliability ®gures are inspected in the modelling data by diameter (Fig. 1a1) and height (Fig. 1b1) classes, it can be seen that the height±diameter prediction model gives almost unbiased estimates for trees of all sizes. However, the model tends to overestimate the heights of the largest trees. The height predictions are slightly better for the test data than for the modelling data, i.e. the means of the residuals are
Table 2 Statistics for the height±diameter model (Eq. (4))a Parameter c0 c1 c2 c3 c4 c5 s2us4 s2um4 s2ut4 s2e0
Estimate 0.038431 1.014449 0.228008 0.029235 0.120584 0.294715 0.000726 0.001853 0.005930 0.005109
Standard error 0.004691 0.001664 0.040968 0.008730 0.007219 0.023196 0.000149 0.000121 0.000187 0.000059
a 2 sus4 is the variance estimate for random stand effects, s2um4 the variance estimate for random measurement effects, s2ut4 the variance estimate for random tree effects, and s2e0 the variance estimate for random errors of the model.
K. EerikaÈinen / Forest Ecology and Management 175 (2003) 355±366
361
Fig. 1. Residual means (-^-, Eq. (7)) and deviations ( , Eq. (9)), and estimates of standard errors of residual means (- - -, Eq. (11)) in the modelling (a1, b1, c1) and test (a2, b2, c2) data with a respect to tree diameter at breast height in 3 cm diameter classes (a1 and a2), estimated height in 2 m height classes (b1 and b2), and stand age in 2-year age classes (c1 and c2). Residuals are obtained for the anti-logarithmic, i.e. back-transformed and bias corrected form of the height±diameter model (Eq. (6)).
smaller, with the exception of the smallest trees (Fig. 1a2 and b2). The means and deviations of residuals analysed against the stand age in the modelling data reveal that the height±diameter model gives unbiased height estimates for trees up to 22 years in age but overestimates the heights beyond 22 years (Fig. 1c1). The height±diameter model also gives slight overpredictions for the trees of the oldest stands in the test data (Fig. 1c2).
Further, reliability ®gures (Eqs. (7)±(10)) were calculated for all observations in the Zambian modelling and test data and also for the Zimbabwean modelling data. Absolute and relative biases for the Zambian modelling and test data were zero. The RMSEs were 1.53 m (8.6%) for the modelling data and 1.42 m (7.5%) for the test data. For the Zimbabwean modelling data the bias was 0.1 m ( 0.4%) and the RMSE was 1.19 m (10.4%).
362
K. EerikaÈinen / Forest Ecology and Management 175 (2003) 355±366
Fig. 2. Relationship between the predicted height±diameter curve (Ð), which is obtained using Eq. (6), and measured dominant height± diameter points (&) (a1 and b1), and the predicted height±diameter curve (Ð) and observed height±diameter points (a2 and b2). a1 and a2 are obtained for the PSP from the tree plantation of Lamba in Zambia; the observed tree diameters and heights are measured at the age of 7.7 (^), 13.8 ( ) and 30.1 () years, respectively. b1 and b2 are obtained for the CCT-treatment from the tree plantation of Martin in Zimbabwe; the observed tree diameters and heights are measured at the ages of 2.2 (^), 4.2 ( ), 8.5 (), 12.4 (^), 18.4 (&) and 25.8 (~) years, respectively.
The development of the predicted height±diameter pattern was illustrated with a PSP-trial from Zambia (Fig. 2a1 and a2) and a CCT-treatment from Zimbabwe (Fig. 2b1 and b2). Fig. 2 shows that both the temporal form and asymptotic development of the height±diameter curve is attained with the model presented. The predicted height±diameter curves are steeper for the earliest stages of stands when compared to the more
even curves obtained for the later measurement occasions. 5. Discussion The aim of this study was to develop a prediction model for the height±diameter pattern of P. kesiya to
K. EerikaÈinen / Forest Ecology and Management 175 (2003) 355±366
be applied in yield simulators and basic calculations of forest inventories. Thus the model predictors used are variables typically measured in forest inventories and also commonly used in yield simulators. Special attention was paid to the modelling of the temporal development of the height±diameter pattern, for which the re-measured sample plot data were very appropriate. According to the results, the presented model predicts both the temporal form and asymptotic development of the height±diameter relationship (see Fig. 2). Steeper height±diameter curves, obtained for younger stands, develop towards more even curves of older stands. In addition, the predicted curves follow the development of stand dominant trees (see Fig. 2). These are the properties required for height± diameter models, which are used in yield simulators based on the predicted development of the stand dominant height and diameter distribution. When the means (Biases) and standard deviations (RMSEs) of model residuals were inspected in diameter (Fig. 1a1 and a2), height (Fig. 1b1 and b2) and age (Fig. 1c1 and c2) classes of the modelling (Fig. 1a1, b1 and c1) and test (Fig. 1a2, b2 and c2) data, it was observed that, even if the variations of classi®ed residuals (RMSEs) are relatively high, the model gives almost unbiased estimates for trees of all sizes. However, the height±diameter model tends to overestimate the heights of the largest trees and gives slight overpredictions for trees of the oldest stands. The overpredictions in the modelling and test data are most probably caused by an insuf®cient number of observations from old stands. In most studies on growth modelling the temporal tree and stand development is linked with the development of dominant height (e.g., Borders, 1989; SaramaÈki, 1992; Hynynen, 1995). This is justi®ed because the dominant height is a measurable stand characteristic and it indicates the site quality in terms of the stand growth rate and yield capacity. It can therefore be used as a predictor of other tree and stand variables needed for the stand growth and yield predictions (e.g., Clutter et al., 1983). With most tree species there is a strong allometric relationship between tree diameter, height and stem volume (e.g., Crow and Schlaegel, 1988). Therefore, prediction errors of stem volumes can be decreased signi®cantly by using the diameter±height ratio as the determinant of the volume (e.g., Laasasenaho, 1982).
363
In fact, the prediction of tree heights is an essential stage when tree taper models are used to obtain stand level volume yields of different timber assortments. However, it is important to consider the risk of using an explanatory variable that contain prediction error, which is obvious in the applications of standard stem volume equations when predicted tree heights are used as explanatory variables of volume models (see e.g., Gregoire and Williams, 1992). There are some earlier studies in which the modelling of the height±diameter relationship has been based on random parameter models, i.e. mixed models, with random and ®xed parameters (Lappi, 1997; HoÈkkaÈ, 1997; Miina et al., 1999). In the estimation of tree height curves, the development of the height± diameter relationship and its linking with the growth of dominant trees is easier to obtain using conditioned forms of nonlinear models (Hynynen, 1995; Miina et al., 1999). Thus, if linear regression estimators and model types are not required, it is possible to use more complex nonlinear model forms and conditions for parameters in the modelling. Parameters are then estimated using nonlinear regression estimators. However, the main aim of this study was to analyse and de®ne a height±diameter model that can be estimated using linear estimators and used in growth and yield simulators based on the development of dominant trees. If the data do not contain spatial or temporal correlation structures, i.e. observations are independent and noncorrelated, the ordinary least squares (OLS) can be used in the estimation of model parameters. However, if these assumptions are false then it is more recommendable to use GLS as an estimator (e.g., Lappi, 1986, 1991, 1997). Nevertheless, the selection of the regression estimator depends on the model forms chosen and, especially, assumptions about error terms and model structures. According to the signi®cance of the variance estimates for the random stand and tree effects of the height±diameter model, it was very evident that observations were correlated spatially. There were also temporal correlation structures caused by re-measurements. Because the basic assumption about noncorrelated residuals did not hold, the random parameter modelling approach with the linear type of model was applied (e.g., Lappi, 1991, 1997; Searle, 1987; Goldstein, 1995; Goldstein et al., 1998). Because of the repeated measurements of
364
K. EerikaÈinen / Forest Ecology and Management 175 (2003) 355±366
the sample plots the hierarchical data structure was broken down, resulting in the cross-classi®ed structure of the data. However, the data enabled the estimation of variances for the random stand, time and tree effects. The height±diameter model of this study is dependent on the dominant height±diameter relationship and its development. The dominant height±diameter relationship is used to scale the curves that are obtained with the model but these curves do not necessarily pass exactly through the observed dominant height±diameter point even if they approximate it, which is logical. Further, the ®xed part of the random parameter model (Eq. (4)) can be used when measured or predicted independent model variables are available. If the model presented is used, e.g. in forest inventories, and if assessed heights of sample trees and stand characteristics, i.e. independent model variables, are available, the model can be calibrated by applying the standard linear prediction theory. Then parameters for random effects are obtained using the best linear unbiased predictor (BLUP) (e.g., Searle, 1971, 1987; Henderson, 1975; Lappi, 1986, 1991; Lappi and Bailey, 1988). Special attention was paid to the selection of the independent variables of the height±diameter prediction model. In this study, the stand age, and dominant height and diameter were used to guarantee the logical temporal development of the height±diameter pattern. It is obvious that the parameters of Eq. (4) are highly correlated. Thus even small changes in parameterisations, i.e. transformations and combinations of predictors, were re¯ected directly in the signs and signi®cance of other parameters estimated. Finally, a simple model form was sought aiming at logical signs and signi®cant estimates of the eventual model parameters. In many studies the logical behaviour of the model with small trees is guaranteed by adding a constant 1.3 to the right hand side of the equation (e.g., NaÈslund, 1937; HoÈkkaÈ, 1997; Forss et al., 1998). By doing this, the model attains height values that are equal to or greater than 1.3 m in every circumstance. In the case of height±diameter equations, which can be transformed into a linear form by taking a logarithm, the constant value, i.e. 1.3, is transferred to the left hand side of equation before taking the logarithm. According
to the observations of Lappi (1991, 1997), however, ln h 1:3 is a very unstable dependent variable for small trees, which degreases estimation properties of a height±diameter equation when data contains measurements of small trees. This was also evident in this study. When the model presented (Eq. (6)) was analysed using the data of this study, it was observed to perform fairly enough with small trees even without the above standard conditioning of the model. The lowest height prediction 1.0 m was obtained for the tree with a diameter at breast height of 0.8 cm. In general, predicted heights attain values close to 1.3 m when the diameter at breast height approaches 1 cm. Finally, the moment when most of trees attain the height of 1.3 m is at the early stage of the rotation, which is a less important period to be considered in stand growth and yield predictions. In the reliability analysis of this study it was observed that the power type of height±diameter model (Eq. (6)) is a reliable and stable predictor. The model has a basic form that smoothes the temporal development of the height patterns of P. kesiya. In addition, the site effect is embedded in the model via the dominant height±diameter dependency of the presented random parameter height±diameter model. These results lead to the recommendation that the model of this study be used as a predictor of the height±diameter relationship of managed stands of P. kesiya plantations in Zambia and Zimbabwe. The height±diameter model presented is an applicable tool with forest inventories and, in addition, with the growth and yield simulators based, e.g., on the predicted development of diameter distribution and dominant height. Acknowledgements I wish to thank Professor Timo Pukkala, Dr. Juha Lappi and Dr. Jari Miina for their valuable comments on this study. The ®eld assessments and data collection in Africa was ®nanced by the Finnish Society of Forest Science and the University of Joensuu, Faculty of Forestry. Financial support by the Finnish Foresters' Foundation is also gratefully acknowledged. The following persons and institutions have supported me in the data collection and ®eld assessments: Mr. Progress M. Sekeli (Division of Forest Research, Zambia),
K. EerikaÈinen / Forest Ecology and Management 175 (2003) 355±366
Mr. Danaza Mabvurira (the Forestry Commission of Zimbabwe), Dr. Ladislaus Nshubemuki (the Tanzania Forestry Research Institute) and Dr. Jussi SaramaÈki. This study is a part of the INCO DC cooperation project Tree Seedling Production and Management of Plantation Forests. Finally, I wish to thank Dr. Lisa Lena Opas-HaÈnninen for revising my English. References Baskerville, G.L., 1972. Use of logarithmic regression in the estimation of plant biomass. Can. J. For. Res. 2, 49±53. Borders, B.E., 1989. System of equations in forest stand modelling. For. Sci. 35, 548±556. Clutter, J.L., Fortson, J.C., Pienaar, L.V., Brister, G.H., Bailey, R.L., 1983. Timber Management: A Quantitative Approach. Wiley, New York, 333 pp. Control Plan 5/3/1, 1969. Systems of measurements in commercial compartments. A System of Continuous Inventory for Pinus Species Established by Industrial Plantations. Mimeograph, Zambia Forest Department, Division of Forest Research, Zambia, 5 pp. Crow, T.R., Schlaegel, B.E., 1988. A guide to using regression equations for estimating tree biomass. North. J. Appl. For. 5, 15±22. Flewelling, J.W., de Jong, R., 1994. Considerations in simultaneous curve ®tting for repeated height±diameter measurements. Can. J. For. Res. 24, 1408±1414. Flewelling, J.W., Pienaar, L.V., 1981. Multiplicative regression with lognormal errors. For. Sci. 27, 281±289. Forss, E., Maltamo, M., SaramaÈki, J., 1998. Static stand and tree characteristics models for Acacia mangium plantations. J. Trop. For. Sci. 10, 318±336. Goldstein, H., 1986. Multilevel mixed linear model analysis using iterative generalized least squares. Biometrika 73, 43±56. Goldstein, H., 1989. Restricted unbiased iterative generalized leastsquares estimation. Biometrika 76, 622±623. Goldstein, H., 1995. Multilevel Statistical Models, 2nd Edition. Institute of Education, University of London, 178 pp. Goldstein, H., Rasbash, J., Plewis, I., Draper, D., Browne, W., Yang, M., Woodhouse, G., Healy, M., 1998. A User's Guide to MLwiN, Multilevel Models Project. Institute of Education, University of London, 140 pp. Gregoire, T.G., Williams, M., 1992. Identifying and evaluating the components of non-measurement error in the application of standard volume equations. Statistician 41, 509±518. Hasenauer, H., Monserud, R.A., Gregoire, T.G., 1998. Using simultaneous regression techniques with individual-tree growth models. For. Sci. 44, 87±95. Henderson, C.R., 1975. Best linear unbiased estimation and prediction under a selection model. Biometrics 31, 423±447. HoÈkkaÈ, H., 1997. Height±diameter curves with random intercepts and slopes for trees growing on drained peatlands. For. Ecol. Manage. 97, 63±72.
365
Huang, S., Titus, S.J., 1999. An individual tree height increment model for mixed white spruce-aspen stands in Alberta. Can. For. Ecol. Manage. 123, 41±53. Hynynen, J., 1995. Predicting the growth response to thinning for Scots pine stands using individual-tree growth models. Silva Fenn. 29, 225±246. Kilkki, P., Siitonen, M., 1975. Simulation of arti®cial stands and derivation of growing stock models from this material. Acta For. Fenn. 145, 33, (in Finnish with English summary). Laasasenaho, J., 1982. Taper curve and volume functions for pine, spruce and birch. Commun. Inst. For. Fenn. 108, 74. Lappi, J., 1986. Mixed linear models for analyzing and predicting stem form variation of Scots pine. Commun. Inst. For. Fenn. 134, 69. Lappi, J., 1991. Calibration of height and volume equations with random parameters. For. Sci. 37, 781±801. Lappi, J., 1993. Biometrical Methods for Forest Sciences. Faculty of Forestry, University of Joensuu. University Press, Joensuu, Finland. Silva Carelica 24, 182 pp. (in Finnish). Lappi, J., 1997. A longitudinal analysis of height/diameter curves. For. Sci. 43, 555±570. Lappi, J., Bailey, R.L., 1988. A height prediction model with random stand and tree parameters: an alternative to traditional site index methods. For. Sci. 34, 907±927. Lappi, J., Malinen, J., 1994. Random-parameter height/age models when stand parameters and stand age are correlated. For. Sci. 40, 715±731. Loetsch, F., ZoÈhrer, F., Haller, K.E., 1973. Forest Inventory: 2. BLV Verlagsgesellschaft, MuÈnchen, Germany, 469 pp. Matney, T.G., Sullivan, A.D., 1982. Compatible stand and stock tables for thinned and unthinned loblolly pine stands. For. Sci. 28, 161±171. Miina, J., Maltamo, M., EerikaÈinen, K., 1999. Modelling the growth and yield of Pinus kesiya in Chati, Zambia. In: Pukkala, T., EerikaÈinen, K. (Eds.), Growth and Yield Modelling of Tree Plantations in South and East Africa. Research Notes 97. Faculty of Forestry, University of Joensuu, Finland, pp. 167±181. Mushove, P.T., 1992. Thirty years of correlated curve trend management trials. In: Piearce, G.D., Shaw, P. (Eds.), Proceedings of the Anniversary Seminar on Forestry Research in Zimbabwe, Mutare, The Forestry Commission, Harare, Zimbabwe, August 27±31, 1990, pp. 75±83. NaÈslund, M., 1937. SkogsfoÈrsoÈksanstaltens gallringsfoÈrsoÈk i tallskog. Rep. Swed. Inst. Exp. For. 29, 1±169, (in Swedish with German summary). Nelder, J.A., 1962. New kinds of systematic designs for spacing experiments. Biometrics 18, 283±307. NyyssoÈnen, A., MielikaÈinen, K., 1978. Estimation of stand increment. Acta For. Fenn. 163, 40, (in Finnish with English summary). Pukkala, T., Pohjonen, V., 1993. Yield of Cupressus lusitanica in Ethiopia. Silva Fenn. 27, 195±207. Pukkala, T., SaramaÈki, J., Mubita, O., 1990. Management planning system for tree plantations. A case study for Pinus kesiya in Zambia. Silva Fenn. 24, 171±180. Pukkala, T., Miina, J., Kurttila, M., KolstroÈm, T., 1998. A spatial yield model for optimizing the thinning regime of mixed stands
366
K. EerikaÈinen / Forest Ecology and Management 175 (2003) 355±366
of Pinus sylvestris and Picea abies. Scand. J. For. Res. 13, 31± 42. SaramaÈki, J., 1992. A growth and yield prediction model of Pinus kesiya (Royle ex Gordon) in Zambia. Acta For. Fenn. 230, 68. Schumacher, F.X., 1939. A new growth curve and its application to timber-yield studies. J. For. 37, 819±820.
Searle, S.R., 1971. Linear Models. Wiley, New York, 532 pp. Searle, S.R., 1987. Linear Models for Unbalanced Data. Wiley, New York, 536 pp. Stage, A.R., 1975. Prediction of height increment for models of forest growth. USDA For. Serv. Res. Paper INT-164, 20 pp. van Laar, A., AkcËa, A., 1997. Forest Mensuration. Cuvillier Verlag, GoÈttingen, 418 pp.