Distance-dependent Models for Predicting the Development of Mixed Coniferous Forests in Finland

Distance-dependent growth models and crown models, based on extensive material, were built for Scots pine and Norway spruce growing in a mixed forest. The crown ratio was also used as a predictor in a diameter growth model to better describe the thinning reaction. The effect of crown ratio on the growth dynamics was studied in simulation examples. Monte Carlo simulation was used to correct the bias caused by nonlinear transformations of predictors and response. After thinnings the crown ratio as a predictor was found to be a clear growth-retarding factor. The growth retarding effect was stronger among pines with thinnings from below, whereas the estimated yield of spruces over rotation was slightly greater when the crown ratio was included than without it. With each type of thinning the effect of crown ratio on pine growth was almost the same, but the growth of spruces was clearly delayed when the stand was thinned from above. Simulation examples also showed that it is profitable to raise the proportion of spruces during rotation, since spruces maintain the growth more vigorous at older ages. The total yield during 90 years rotation was about 20 % higher if the stand was transformed into a pure spruce stand instead of pine.


Introduction
One of the most important aims of forest management is to increase economic return.To maximize utility, the forest owner frequently maximizes the monetary income.By integrating individual forest owners' utilities, maximal total utility of the forest management can be reached.This strong supposition assumes e.g. that prices of timber assortments and interest rate are determined correctly.Maximum economic return can be reached by optimizing the treatment of the single stands.Because numerous different treatment schedules cannot be compared in real stands, simulation seems to be a reasonable alternative to study the effects of various treatments.The optimal treatment denotes proper tree species composition in the initial stand, correct timing of cuttings and correct tree selection in cuttings.The more accurately development of the stand can be predicted, the more reliable are the optimum treatment programmes obtained in simulations.Therefore, by increasing the accuracy of simulation models, the reliability of stand development predictions can be increased; and along with that, the more profitable forest management decisions can be put into practise.
In a few studies it has been shown that on sites of medium fertility it is profitable to grow mixed conifer forests (e.g.Jonsson 1962, Kerr et al. 1992and Pukkala et al. 1998), and in certain cases thinning from above is more advantageous than thinning from below (e.g.Vuokila 1970and 1977, Hynynen and Kukkola 1989, Mielikäinen and Valkonen 1991and Vettenranta 1996).When thinning from above is simulated, some problems arise.Competition affects small trees more than large trees because a larger tree can use limited growing resources more efficiently than a smaller tree; therefore, the change in competition after thinning will be emphasized if large trees are removed instead of small ones.Thus, predictions of growth will be less accurate.
Several methods can be used to increase the accuracy of growth predictions, when either the structure of initial stand or applied treatments deviate from the common practise.Distance-dependent growth models are generally used to describe the effect of competition on growth of a single tree.Especially in mixed forests, competi-tion indices are useful tools for describing the variation in stand structure so that the development of single trees can be predicted more precisely (e.g.Daniels et al. 1986, Pukkala 1989, Pukkala and Kolström 1991, Biging and Doppertin 1992, 1995).However, some problems arise when the spatial structure of a stand suddenly changes in connection with thinning.
Stand basal area and competition indices are changeable external growth factors that can be used as predictors in tree-level growth models.Real growth probably does not increase as rapidly as growth prediction does in a single-tree growth model, when external growth factors are changed.Real growth depends partly on the needle mass of the tree, which does not change as rapidly as competition after thinnings.Therefore thinning reaction can be made more natural by including the needle mass in the growth model, e.g. with the assistance of crown ratio (height of the living crown divided by total height of the tree).
Crown size has proved to be a useful indicator of thinning reaction (e.g.Smith 1986 andShort andBurkhart 1992).In Finnish conditions the effect of crown size has been studied by e.g.Hynynen (1995a and1995b) in pure Scots pine (Pinus sylvestris L.) stands and by Mielikäinen (1985) in mixed stands of Norway spruce (Picea abies (L.) Karst.) and birch (Betula pendula Roth.and Betula pubescens Ehrh.), in which the crown ratio significantly improved the accuracy of the models.According to the study of Hynynen (1995b), the model with crown ratio as independent variable resulted in underestimation of the diameter growth in the first 5-year period after thinning, when stands were heavily thinned (more than 50 % of the basal area removed); but in unthinned and moderately thinned stands the predictions were satisfactory.
In general, the crown ratio is predicted either by an allometric (e.g.Dyer and Burkhart 1987) or an incremental model (e.g.Short and Burkhart 1992).According to the study of Liu et al. (1995), the methods do not differ markedly in their ability to predict the crown size.However, allometric models can be utilized, although the simulation input does not include information on crown size.Both non-linear (e.g.Hynynen 1995a andLiu et al. 1995) and linear (e.g.Mielikäinen 1985) ap-proach have been used in modelling crown ratio, expressed as a variable between 0 and 1.
The aims of this study were (1) to build distance-dependent models for development of mixed conifer forests in order to analyze the effect of crown ratio on growth predictions and (2) to verify that the development of a mixed conifer stand can be predicted reliably over the rotation.

Material
The study material consisted of about 350 study plots measured by The Finnish Forest Research Institute throughout Finland, except for the two northernmost districts of the Board of Forestry (Fig. 1).The study plots were permanent study plots on mineral soil, located systematically according to the sample areas of the National Forest Inventory: in each sample area one stand was chosen randomly.Stands with notable damage were not entered to the data.The data on sapling forests, dominant height lower than five meters (TINKA-data), were measured two times during 1984-1989 and the data on older forests, dominant height over five meters (INKA-data), were measured three times during the years 1976-1992.The intervals between measurements were five years (Gustavsen et al. 1988).
Most of the stands are of natural origin.About 35 % of the trees in the TINKA-data have been planted and about 15 % have been seeded.In the INKA-data the origin was not included and it is probable that most of trees are of natural origin.
In each stand, three circular study plots were established including at least a total of 100 trees.The distance between the study plots was 40 meters.The area of each stand compartments was at least 0.5 hectare.The areas of the study plots varied between 0.002 and 0.023 ha in the TINKA-data and between 0.008 and 0.13 ha in the INKA-data.The mean areas were 0.006 and 0.011 ha, respectively.Sample plots (area one third of the study plot) were located and centered within each study plot.A buffer zone, the difference between a study plot and a sample plot radius, varied between 0.26 and 8.46 m, with a mean of 3.68 m.All stand compartments were measured for distance from the equator and mean age at breast height.The mean age was not determined as the mean age weighted by the basal area, but as the average age of the dominant trees at breast height (increment cores or an estimate).The time from the last thinning was also observed.Forest site types were determined according to the instructions of the National Forest Inventory (Valtakunnan...1977).From the total INKA-and TINKA-data were selected those coniferous mixtures in which the quantity of mixed species was at least 10 % of the number of stems (Fig. 2).After that, the TINKA-data consisted of about 16 000 and the INKA-data about 120 000 observations.About 5 % of the observations in the diameter growth modelling were birches (Betula sp.) and they were treated as pines when competition indices were calculated.
All trees on the study plots were measured for species, location and diameter at breast height (referred to here as diameter).Trees on sample plots were also measured for height, and some sample trees were measured to determine the height of the lower edge of the crown.Lonely living branches were not considered as a part of the crown if two or more dead whorls of branches were found above them.The plot measurements were used to prepare the model for tree height.Näslund's (1936) curve (h = 1.3 + d 2 /(a + bd) 2 , in which h is tree height (m) and d is diameter) was fitted to the data: thus the height of the non-sample trees could be predicted.The dominant height of a stand was determined as the mean height of 100 of the thickest sample trees per hectare.Only those trees for which differential growth, height growth and crown height were measured were used as data sources for the estimated model in question.The differential growth in diameter and height was combined with earlier measurements, so the future five-year growth could be used in modelling.
The maximum number of observations in modelling was about 10 000 pines and 11 000 spruces, and the number of trees was slightly more than half the number of observations.About half of the pines were located in a forest with fresh mineral soil or on better sites (e.g.Myrtillus site type (Cajander 1909)), while 87 % of the spruces were situated on sites that were equally good.Spruces were located in stands that were, on average, older and more heavily stocked than the pine stands were.Spruces were also more numerous in the southern and western parts of the country (Table 1).

Modelling
The ordinary least-squares (OLS) method for fitting the regression model presupposes that all observations used in modelling are independent.
The data in forest research do not always meet this requirement.If, for example, several trees are measured in the same stand or there are many observations for the same tree at different times, this intra-unit correlation and hierarcical data structure can be taken into consideration by using the generalized least-squares (GLS) method for fitting the model (Goldstein 1995).Better estimates than OLS estimates can be obtained by GLS estimation if the covariance structure of the data is considered.The covariance matrix can be estimated iteratively by starting with OLS estimates for the parameters of the function in ques-tion.The residual variation can be divided into components described by random parameters (Lappi 1986, Lappi andBailey 1988).The random parameters are assumed to have zero mean, fixed variance and they are assumed to be uncorrelated with each other (Searle 1987).However, tree dimensions within a stand can be spatially correlated (Penttinen et al. 1992), due to, e.g.competition between trees.Thus competition indices between nearby trees also are correlated.GLS estimation, in the way it has been used in this study, does not take this correlation into consideration.However, it was expected to give better estimates than OLS estimation, because the correlation between competition indices was assumed to be relatively small.In this study three random variables were used: stand variable (observations in the same stand are correlated), tree variable (observations from the same tree at different times are correlated) and error variable.
So the model will be in the following form: in which y ijk is 5-year diameter or height growth or crown ratio of tree i in plot j in measurement time k, X1 j ,...,Xm j are stand characteristics of plot j, which are not expected to change when the moment of time changes, X(m + 1) jk ,...,Xn jk are stand characteristics of plot j, which change during the measurement time k, C(p) ijk ,C(s) ijk are competition indices computed for pines and spruces, respectively; d ijk , h ijk and cr ijk are diameter, height and crown ratio of tree i in plot j and measurement time k, respectively; p j is a random plot factor, t ij is a random between-tree factor and e ijk is a random within-tree factor.

Competition Indices
The structure and the growth dynamics of mixed forests are more complex than those of pure stands, due to different ecological properties of species (growth rhythm, demand of light, use of nutrients etc.).Therefore, a spatial arrangement of species and individual trees must be known so that between-trees interaction can be taken into consideration when single-tree growth models are used for the detailed description of a stand structure development.The knowledge of spatial distribution of trees is crucial specially when the structure of a stand or applied treatments are exceptional.One way to take the spatial arrangement into consideration is to use competition indices to describe the interaction between trees.
Competition indices based on vertical angles defined by a horizontal plane to the tree top inside a certain radius (Pukkala et al. 1994) have proved quite useful in coniferous forests, although some problems arise.There is no implicit (natural) reason why competition would suddenly stop at a certain distance, but it is more realistic to assume that the competition decreases as the distance from a subject tree increases.Secondly, in mixed coniferous forests the dominated tree storey often consists of spruce, due to its better adaptation to low light conditions.If there is a marked difference between tree heights, the horizontal angle does not increase with the competition.Then only the number of competitors affects the competition.
In this study several competition indices, based on the same type of idea, were tested.First, there was no specified limit to the distance at which the competition would stop, but this was determined by the stand structure.Secondly, the smaller the trees were, the more quickly the effect of competition decreased along the distance, for example, due to the width of the root systems.Thirdly, competition indices had to be appropriate for the mortality model.The best applied competition indices, in respect to their ability to explain the residual variation with other applied variables, were the sum of the proportion of heights of those competitors exceeding the paraboloid surface centered on the subject tree from a certain relative height (Tables 2 and 3, Fig. 3), as follows: in which n is number of trees; h j is height of the competitor j; h i is height of the subject tree; D ij is 57 distance from competitor j; H d is dominant height of the stand; a is 0.5, 1 or 2; and b is 0.4 or 0.8

Mortality Models
In order to simulate forest stand dynamics, a mortality model is usually needed.Mortality models are often based on stand variables and the relative size of a subject tree (e.g.Ojansuu et al. 1991).That kind of model does not consider the spatial distribution of trees.In particular, thinning from above may change the spatial distribution of trees, making it unfavorable.In this study the mortality models were based on the competition indices so that if the competition compared to height of a tree exceeded a certain level, the tree was removed as a mortality.This approach allows that in regular stands the growing-stock density may become higher than that in clustered stands.

Diameter Growth Models
The range of variation of the tree size along with the stand age was so broad in the material used to prepare the models that most of the independent variables had to be modified to linearize the dependence.The models were not only modified with the assistance of statistical criteria, but also by studying how diameter and height growth models function together.The relative size of a subject tree (or some modification) is one of the most significant predictors in each model.A certain modification was chosen not only because of the ability to predict the growth, but also according to the stem form (d/h) development during the rotation.Therefore predictors are not necessarily the best ones from a statistical point of view.The square-root transformation was used in most models to decrease the heteroscedastici- and the following for spruce, The following diameter growth models were estimated for pine, in which id is future 5-year diameter growth (mm); d is diameter (mm), h is height (m); C p 8, C s 6 are competition indices calculated from pine competitors and from spruce, respectively (see Formula 1 and Table 2); Nc is distance from the equator (km); Tg is stand age at breast height (a); G is stand basal area (m 2 ha -1 ); dummy4 is 1, if the stand was located on a worse than the fresh mineral soil forest site (e.g.Vaccinium site type), otherwise 0; p is a random stand factor; t is a random between-tree factor and e is a random within-tree factor.Subscript i refers to a tree, j to a plot and k to a time of measurement.The tvalue of parameter estimates (the parameter estimate divided by the standard error of the estimate) and the z-value of random parameters (the variance estimate divided by the standard error of the variance estimate) are presented in brackets below the parameters.
The within-tree residual mean square error was estimated to 0.427 for the pine model and 0.318 for the spruce model.The corresponding be-tween-stand and between-tree components of the total residual variance were 0.124 and 0.189, respectively, for the pine model and 0.226 and 0.214, respectively, for the spruce model.The estimated pure error at tree level (e ijk ) accounted for about 57 % of the total residual variation in the pine model and 42 % in the spruce model.The degree of determination of the mixed model part was 0.58 for the pine model and 0.57 for the spruce model.Relative RMSE's (RMSE divided by mean value of the predicted variable, Eq. 5) were 0.223 and 0.249, respectively (Table 4).

Diameter Growth Models with Crown Ratio
Diameter growth models with the crown ratio were as follows for pine, in which cr is crown ratio (height of the living crown divided by total height of the tree).Other predictors are as in Equations 3 and 4. By including the crown ratio in the models, the total residual mean square error decreased about 10 % for each model.The between-tree component of the total residual variance of the pine model decreased about 35 %.In the pine model the crown ratio is the most significant predictor, in the spruce model it is the second most significant.In the pine model the competition of spruce and stand age alone are no more significant predictors.The number of observation decreased from 8 900 to 4600 in the pine model and from 11 000 to 4200 in the spruce model, when crown ratio was included.
Both the growth model without and with the crown model gives the same results if the crown ratio is at a certain level.That certain crown ratio decreases with the age of the stand and is higher with the model of spruce than with the model of pines (Fig. 4).
The diameter growth models (Eqs.3, 4 and 6, 7) give the same results when the crown ratio is 1 for spruces and slightly less than 0.7 for pines at the stand age of 5 years at breast height.At the stand age of 80 years, the crown ratios are about 0.75 and 0.45, respectively.When the crown ratio was decreased about 10 percentage units, the diameter growth also decreased about 10 %, and vice versa.4. Crown ratio when diameter growth models without and with crown ratios gives an the same diameter growth (solid line) and crown ratio when diameter growth with crown ratios gives 10 % smaller and greater diameter growth (dashed line).
If crown ratio is greater the diameter growth is greater and vica versa.Spatially the stand is regularly distributed and all trees have the same mean dimensions according to the stand age according to Koivisto 1959.

Models for Crown Ratio
In order for crown ratios to be predicted for the diameter growth model, static crown ratio models were built for pine, development of pines, as in the diameter growth models.On the other hand, in both models the stand basal area is one of the most significant predictors, and in the pine model the competition described by the stand basal area affects the crown ratio more than the competition of each species separately (Fig. 5).

Height Growth Models
Height growth models were as follows for pine, The crown ratio for pine can be predicted much better than that for spruce.The relative RMSE's were 0.17 in the pine model and 0.33 in the spruce model.In both models the stand basal area and the competition are the only factors that decrease the crown ratio.The competition of pine is not a significant predictor for the crown ratio of spruces, but for pines the competition of each species decreases the crown ratio.This aspect indicates that pine is not as harmful a competitor for spruces as spruce itself, but the competition of each species significantly affects the in which ih = future 5-year height growth (cm).
The between-tree component of the total residual variance was 5 % in the spruce model, while in the pine model it was about 70 %, because it was not possible to determine the within-tree component; but the between-tree component includes also the within-tree variation.Pines reach maximum height growth at a height of 4 meters, while the height growth of spruces reaches its maximum at a height of 5 meters.According to Koivisto (1959), the biological age of a pine stand is then about 20 years and the age of a spruce stand is about 27 years.According to the competition indices only competition by its own species decreases the height growth of each species; on the other hand, in the pine model the difference between tree height and dominant height promotes the height growth.The degree of determination of the fixed part of the model was 0.63 for the pine model and 0.37 for the spruce model, even though the relative RMSE's were 0.37 and 0.27, respectively.

Mortality Models
The idea of the crown model is to describe the vitality of a tree.Unfortunately, the crown ratio models are not flexible enough to be used as mortality models.The height or diameter growth, calculated with their respective models, breaks off before the crown model diminishes the crown of a tree minor.The tree crown withers away rather than diminishes little by little to non-existence.As growth models built for pure stands do not give reliable growth predictions in mixed forests, it was probable that existing mortality models would not give reliable predictions in mixed forests.Therefore, new mortality models were built in order to remove single trees under heavy competition and thus avoid unrealistically high stand density.The mortality was determined with the aid of envelope in the present data.Certain competition indices were compared with the difference between a tree height and a dominant height and a limit was observed.It was assumed that trees which have exceeded the limit were dead, and therefore no observation was made.The total number of observations was about 18 000 in pine models and abour 21 000 in spruce models.The mortality models were formed as functions of the limits as follows:

Model Validation
Diameter growth models were tested with independent data sets collected in mixed coniferous forests in North Karelia.One data set consisted of a one-storey conifer mixture (Pukkala et al. 1994) and another consisted of a two-storey mixture (Pukkala et al. 1998) in which pines were the dominant trees and spruces the under-storey.All diameter growth models functioned quite well in the one-storey data set, except that the growth model for pines with the crown ratio (crown ratio in test-data were estimated) gave an average 0.3 cm underestimation of the 5-year diameter growth and the high predictions in each spruce growth model were overestimated.Conversely, in two-storey stands the small predictions of pine models were underestimated and the high predictions were overestimated.The predictions of the spruce models were slightly underestimated.All models were also tested with a partitioned data set.Every fifth observation was separated into a different data set and all models were rebuilt with the remaining 80 % of the observations.The parameter estimates did not significantly change with the subset and the predictions were very good.
The geographical distribution of prediction errors was also studied by viewing the random plot factors, which indicate the growth level of each stand compartment.The diameter growth model for spruce underestimates the growth in Ostrobothnia and in southeastern Finland and overestimates the growth in North Karelia.The pine model overestimates the growth in the area from Nurmes to Vaasa.

Simulation system
To test the function of the growth models and the effect of the crown ratio model, two simulation models were built.Both simulation models contained the same height-growth models (Eqs.10 and 11) and mortality models (Eqs.12 and 13).One model also contained the diametergrowth models without a crown ratio (Eqs.3 and 4), and another simulation model contained diameter growth models with crown ratio (Eqs.6 and 7) and the crown ratio models (Eqs.8 and 9).
The structure of a stand is a crucial factor when the development of a stand is simulated with these models.In order for the growth models to give reliable predictions, the tree dimensions and spatial distribution of the trees within a stand must be correct.Therefore, in the following simulation examples, one real stand was used as the initial stand.This stand was naturally regenerated and nearly even-aged, located in North Karelia, Finland.The area of the study plot was about 0.11 ha and was located on the site of medium fertility (MT, Myrtillus site type (Cajander 1909)).Spruces were, on average, smaller than pines (Table 5).
The simulation followed a pattern.First, a temporary buffer zone was generated around the study plot by assuming that the plot was surrounded by similar plots on all sides.After that competition indices could be calculated.Then the crown ratios were generated with the crown ratio models if they did not exist and they were needed.Diameter and height growths for five years were calculated.The trees with inadequate vitality were removed with the mortality model, which removes trees according to the surrounding competition.Next, the stand characteristics were calculated with the new tree dimensions.The buffer zone was generated and competition indices were calculated again with the new tree dimensions.The new crown ratios were calculated with the new tree dimensions and the stand characteristics.The lower edges of the crowns were determined with the crown ratio.If the earlier lower edge of the crown was higher than the new one, the crown ratio was corrected with the earlier one.The volume of trees and the proportion of timber assortments were calculated with the functions of Laasasenaho (1982).Finally, the buffer zone was removed.
Thinnings were simulated by setting the number of thinnings and a certain basal area limit for each thinning.If this limit was exceeded during the present 5-year growing period, the growth was predicted by interpolating, so that the basal area limit was attained.The trees of each species were divided into three diameter classes by the number of the stems.The different types of thinnings were implemented by giving a wanted proportion of removing trees for each diameter classes.The treeselection method for harvesting is described in detail in the study of Pukkala et al. (1998).
Two different harvesting methods were applied, so that the function of the models could be compared in different situations.For studying the development of the mean diameter and the height of species, the stand was thinned from below, so that 80 % of the trees in the lowest diameter classes of each species were removed, 70 % from the second class, and also 10 % from the highest diameter class were removed to describe to a removal of trees with a heavy competition against them and a bad quality.
The first thinning from below was conducted when the stand basal area reached 30 m 2 ha -1 .The post-thinning basal area was about 20 m 2 ha -1 and the age of the stand was then about 50 years.The second thinning was conducted when the basal area reached 35 m 2 ha -1 and the age of the stand was about 75 years.The post-thinning basal area was about 22 m 2 ha -1 and the final density of the stand was about 400 stems per hectare, while in the initial stand the density was about 1500 stems per hectare.
In most of the examples the thinnings from above were conducted by removing all trees from the highest diameter class of each species.The first thinning was conducted when the stand basal area reached 30 m 2 ha -1 .The post-thinning basal area was about 11 m 2 ha -1 , and the age of the stand was then about 50 years.The second thinning was conducted when the basal area reached 30 m 2 ha - 1 and the age of the stand was about 80 years.The post-thinning basal area was about 14 m 2 ha -1 and the final density of the stand was 550 and 520 stems per hectare depending on the simulation model.The spatial distribution became unfavorable by reason of the thinning from above, and over the rotation the mortality model removed about 160 stems with the model without the crown ratio and about 190 stems with another model.

Correction of Bias
When growth predictions are calculated with non-linear transformated growth models, the predictions can be biased (e.g.Gertner 1991).If the growth is to be predicted over long periods, e.g. over a rotation, so that earlier predictions will be used as the predictors in the growth models, the bias will be accumulated.The longer the simulation period is the more uncertainty will appear in predictions.Therefore, several corrections of bias were used in the simulation models.The bias caused by the transformation of the response and the bias caused by non-linear transformated predictors were corrected by using the Monte Carlo simulation, in which a normally distributed random factor was added to the predictors before the prediction was calculated (Kangas 1996 and1997).This was repeated 1000 times and the final prediction was the mean value of these simulations as follows: and finally in which f(y i' (1) ijk ) is some transformation of the response function at the present time, g is some linear function of the predictors, x i' (1) ijk is predic- tor that is predicted with a model and corrected with the Monte Carlo simulation, e f(yi') is normally distributed random factor whose expectation value is 0 and variance is the variance estimate of the prediction model of the predictor in question (sum of variances of p j , t ij and e ijk in Equation 1), f -1 is inverse function of the transformation, number 1 in brackets refers to the present time and number 0 refers to the previous growing period.Other symbols are as in Equation 1.
Because the structure of the models is hierarchical, all simulated corrections in Equations 14 and 15 can be made at once, so that all the independent random factors are added, the response is calculated, the inverse transformation of the response is made and this was repeated several times and then the final prediction will be calculated as an average of the single predictions (Eq.16).

Development of Dimensions
The development of the dominant height in the mixed forest and in a pure pine forest (treated with thinnings) on the Myrtillus site type, according to Koivisto (1959), are almost the same until the age of 80 years when the large spruces start to determine the development of the dominant height (Fig. 6).In the present simulation models the height growth of spruces does not slow down as soon as in a pure spruce stand; and the growth rhythms of species in the mixed forest are much more alike than in a pure stands.Growth of the mean height of spruces is more vigorous than that of pines, and also the thinnings affect the mean height of spruces more, even though both species were thinned with the same removal percentages.
The mean diameter of each species in the initial stand is slightly bigger than the corresponding mean diameters in pure stands (treated with thinnings) (Koivisto 1959) (Fig. 6).Accordingly, the simulated development of the mean diameter in the mixed forest is slightly more rapid than that in the pure stands.When the simulation model with the crown ratios (SCRY) is used, the development of the mean diameter of spruces is slightly faster in proportion to the diameter growth of pines than in those cases where the simulation model without the crown ratios (SCRN) was used.
The effect of crown ratio on the growth dynamics after thinnings and the division of growth into the diameter classes was studied by applying two different treatment alternatives: the aforementioned thinning from below and a pure selection from above, in which all the trees in the highest diameter class of each species were removed.
Before the first thinning from below, over a period of 15 years, the mean diameter of pines with SCRN had become about 1 mm larger and the mean diameter of spruces about 2 mm smaller than in the other model.After the thinning from below, the difference in mean diameters (Dg without cr -Dg with cr) of pines between the simulation models decreased to 0 mm and the difference in the mean diameters of spruces decreased about 4 mm, being 6 mm after the thinning (Fig. 7).
During the first growing period after the first thinning, the growth of each species was retarded about 4 mm when SCRY was used, as compared to SCRN.In the beginning of the period between thinnings, the crown ratio is clearly a retarding factor.The retarding effect decreases during the time when the crowns become longer, and the growth is more vigorous with SCRY in the end of the period between the thinnings.After the second thinning, the difference between the growths followed the same pattern as after the first thinning.At the age of 90 years, the final difference between the mean diameters obtained with the simulation models was 8 mm with pine and -4 mm with spruce.When the thinnings were performed from above instead of from below, the change in mean diameters of spruce after thinning was the opposite (Fig. 7).After the first thinning, the difference in the mean diameters obtained from simulation models was the same.During the following 15-year period, the crown ratio can be seen as a clear growth-retarding factor for each species.With thinning from above, the effect among pines and spruces is almost the same.
The final difference between the models was 11 mm with pine and 13 mm with spruce for the benefit of SCRN.When different types of thinnings are compared, the effect of a crown ratio on the diameter growth of pines is slightly less when thinning from below; but after the thinning from above for spruce, the retarding effect of the crown ratio is clear.The effect of the crown ratio on growth is obvious after the thinnings regardless of tree species or type of thinning.
The simulated volume yield of the stand is slightly faster than the yield of a pure pine stand on a Myrtillus site type, but slightly slower than the yield of a pure spruce stand on a Oxalis-Myrtillus site type (Koivisto 1959).The volume of the initial stand was clearly lower than the volume of stands in the standard of comparison (Fig. 8).

Effect of Species Composition
The effect of species composition on the volume yield during rotation was studied by adjusting the removal ratio on each of six diameter classes in two thinnings from below and above.The prethinning basal areas were 30 m 2 ha -1 before the first thinning and 35 m 2 ha -1 before the second  thinning.Post-thinning basal areas were about 20 m 2 ha -1 and 22 m 2 ha -1 , respectively.In the following comparisons the rotation was 90 years.
With each simulation model and type of thinning, the total volume yield over the rotation was the greater the greater the proportion of spruce in the end of the rotation (Fig. 9).This indicates the fact that the volume growth of spruces is not delayed as soon as the growth of pine when the stand approaches the final cutting.On the other hand, the volume yield does not increase linearly along with the proportion of spruce, but a mixed forest effect on the volume growth can be seen.
SCRN and the alternative with the thinning from below gave the greatest volume yield with every species composition.In all cases SCRN gave a higher volume yield than SCRY.In the cases in which all the trees are spruce in the end of the rotation, the volume yield is about 100 m 3 greater than in a stand occupied only by pines in the end of the rotation.
In addition to the expected mixed forest effect, mortality causes less productivity in the pinedominant stand.In almost every case, the proportion of pines indicates the number of stems removed by mortality (Fig. 10).Thinning from above also promoted mortality, from 125 to 325 stems ha -1 .By thinning from below, the mortality varied between 75 and 225 stems ha -1 .When SCRY was used, the mortality was lower in all cases with the thinning from below and greater  when thinning was from above.The mortality among spruce after the thinnings from above also indicates that the growth of small spruce suffers from competition.
The function of the mortality models was also studied by growing the stand without thinnings.The density started to decrease from about 1500 stems ha -1 at the age of 45 years when the stand basal area was about 25 m 2 ha -1 (Fig. 11).During the period between 45 and 60 years more spruce were dying and thereafter the mortality among each species was almost the same until the stand age of 100 years.Older pines were still dying, but the number of spruce no longer decreased markedly.The growth of the stand basal area was delayed when it approached the level of 40 m 2 ha -1 at the age of 70 years.The stand basal area remained steady about the level of 40 m 2 ha -1 for 50 years; then at the age of 120 years, it again started to grow slowly.

Model Evaluation
Mixed coniferous forests are very common (about 20 % of the forest area) in Finland (Kuusela and Salminen 1983).Development of a mixed forest cannot be predicted reliably enough with models built for pure stands.This paper presents growth models that allow us to predict the development of the pine-spruce mixture over rotation throughout southern Finland.The form of the models also allows research on the interaction between trees and the effect of thinning methods.The main purpose of these models is to determine whether, by including crown ratio into the growth models, the development of mixed coniferous forests can be predicted precisely enough, in order to determine the optimal treatment schedule in spite of thinning type or structure of the initial stand.
The material used to prepare the models is so comprehensive that it allows in this sense an extensive use of the models: on all natural forest site types of a conifer mixture, fairly varied in terms of tree age and size distributions, from stand age at a breast height of about 5 to over 100 years and anywhere in southern Finland.Development of clearly uneven-aged and twostorey forests cannot be predicted as reliably as that of more uniform forests.Although, the data as a whole does not represent mixed coniferous forests in a stand level, in a single tree level there are observations enough to describe the variation in competition of each species reliably (Pukkala et al. 1994).In the pine data 23 % and in the spruce data 37 % of observations include non-zero competition caused by other species according to the used competition indices.
The sensivitity of the models to inputs was studied by changing the mean value of the predictors one at a time 10% of the range greater and smaller and then by examining the change in the results.Diameter growth was affected most either by the stand basal area (7-8 %) or by the crown ratio (7-11 %), if it was included.The effect of the stand basal area decreased about 25 % in the pine model and about 40 % in the spruce model when the crown ratio was included.Thus the post-thinning diameter growth is retarded in two ways.First, the effect of the suddenly decreasing stand basal area diminishes and secondly, the growth is also retarded by the reduced crown ratio, as it's becoming larger depends on the height growth.
The growth models are affected most by crown ratios and less by diameter (only pine), height, stand age and stand basal area.Although the crown ratios are also affected by the aforementioned variables and they affect in the same direction as in the diameter growth models, except for tree height.The crown ratios are not very sensitive to any variable, and the competition indices have the least effect.When the maximum competition is added, the crown ratio of pine decreases from 56 to 34 % and that of spruce from 78 to 44 %.
Stand age has become emphasized in the spruce height growth model, and it is a very significant predictor in the other models, too.In all the models either the diameter or the height of a tree is related to the stand age as a predictor.So, by determining stand age wrongly compared with tree dimensions in the initial stand, unreliable results will be obtained.
Most stands are of natural origin.About 35 % of trees in the TINKA-data have been planted and about 15 % have been seeded.In the INKAdata the origin was not included and it is probable that most of the trees are of natural origin.So, the possible different development of the cultivated sapling forests is also taken into account to some extent.As long as the structure of a initial stand is determined correctly and there is no special factor that would bias the future development of a stand, there is no particular reason why in most circumstances these models would not give reliable enough predictions to determine the optimal treatment schedules.

Analysis of Results
The results of simulation examples agree mainly with the development of pure pine and spruce stands according to Koivisto (1959).However the development of spruces with the present models follow the development of a pure stand on the better forest site type (Oxalis-Myrtillus site type)(Fig.6).This may be caused by difficulties to determine the forest site type and by the mixedforest effect.The development of the stand basal area also follows mainly the development of a pure spruce culture on the good site class (H 100 = 30) according to Vuokila and Väliaho (1980).But the development is markedly more vigorous than that of a pure pine culture on a corresponding site class (H 100 = 30).
By including the crown ratio into the models, the estimation of the diameter growth is affected noticeably after the thinnings.The thinning reaction (an increase in growth after thinning) of each species is delayed when the crown ratio is included.The retarding effect lasts 10-20 years and after that, before the next thinning, the crown ratio starts to increase the growth compared with the model without the crown ratio, specially when spruces are dealed with thinning from below.By including the crown ratio into the growth models the reliability of growth predictions after thinning from above was expected to increase.However, the data is mainly consisted of stands thinned from below (common practise in Finnish forests).Although, the data includes small single trees which have been released from competition and the (d or h)/Tg ratio and the crown ratio describe the vitality of those trees, the behaviour of the earlier suppressed trees after the thinning from above is still uncertain.Thus, empirical studies about the thinning reaction after thinning from above is needed.
The difference in results between the simulation models (SCRY and SCRN) is expected when the development of the spruce diameter is observed.Thinning from above increases the difference between the models more than thinning from below does, as the crown ratios have become smaller and big spruces in the young stand grow better in the simulations with the crown ratio.The better growth before the first thinning below compensates the delayed development after thinning, so in the end the simulation models produce almost the same mean diameter for spruce.Whereas, by thinning from above the obtained mean diameter is finally clearly greater without the influence of the crown ratio.
The effect of crown ratio on pines is not as clear when the difference between the thinning types is compared.However, the crown ratio decreases the development of mean diameter slightly more when the stand is thinned from above.All aforementioned aspects indicate the fact that crown ratio is a profitable indicator of the thinning reaction and that thinning from above is more suitable treatment method for pine than it is for spruce (e.g.Mielikäinen and Valkonen 1991).By using the simulation model with the crown ratio, the development of the stand after the thinning from above can be predicted more precisely and the results will be more reliable than the predictions without the crown ratio.However, not even the present models will necessarily provide reliable predictions in extreme conditions, because the growth models give predictions of both non-zero height growth and diameter growth, even if the crown ratio would be zero (a dead tree).
The crown ratio of pines that produce the same diameter growth as the model without the crown ratio (Fig. 3) is about the same level as that simulated by Hynynen (1995a) in a thinned pine stand.On average the crown ratio of each species in the present simulations is slightly greater than the crown ratio in the study of Hynynen and Siipilehto (1996), but is smaller than the growth equalizing crown ratio in Fig. 4.
The development of the height structure of the stand in simulations supports the findings of Vettenranta (1994) and Ojansuu et al. (1991), which showed spruces growing much more vigorously than found by Koivisto (1959).Thus, the problem is to determine the yield ability of a certain stand.In Finland the forest site types are traditionally determined with the aid of soil vegetation (Cajander 1909).In mixed forests, problems arise because each species has a different effect on soil vegetation.In pure spruce stands, growth may be delayed by deterioration of the soil.The decomposition of litter is delayed, because of a lower temperature in the ground layer and the acid litter of spruces (Mikola 1954).Different tree species also use nutrients, solar radiation and water differently.The species mixture affects the properties of the humus layer and the nutrient condition of the soil (Sepponen et al. 1979), and along with that the development of the dominant height, which is widely used to describe the yield ability of a stand, is also affected (Tamminen 1993).
When the dominant height-age curves are used in determining a site index, there is also a problem: the development of dominant height depends on the species composition, so the same stand looks like better growing place for spruce than for pine (c.f.Nigh 1995).The problem already arise when an attempt is made to deter-mine the stand age.Different species reach breast height at different ages (Gustavsen 1980), so determination of a site index is uncertain.Even in planted stands the difficulties occur; and the problem is more complicated in naturally regenerated stands where the structure of the previous stand affects the sapling forest.For aforementioned reasons both site classifications in an onedimensional form based on dominant height or on soil vegetation seem to be inappropriate for mixed forests.
The present simulation models were build for use throughout southern Finland.However, the simulation examples indicated that the development of a stand located on a particular latitude cannot be simulated by merely changing the location variable, as the forest site type cannot be described only by a dummy variable, but that the structure of a initial stand is crucial for correct predictions.
Competition indices are good tools for describing the interaction between tree species.According to the competition indices of the present models, spruce suffers very little from the competition with pine.A finding that agrees with that of Pukkala et al. (1994).Pine growth decreases along with the competition of each species, but when crown ratio was included in the diameter growth model, the competition of spruce became not significant.Competition is taken into account with the aid of the crown ratio, which is affected by the competition of each species.On the other hand, the areas of the study plots were so small that the stand basal area often describes the total competition around a subject tree.Competition indices only calibrate the competition when the species composition, the spatial arrangement or the density around the subject tree varies.The effect of the competition indices is emphasized when the spatial arrangement of a stand is unfavorable and some trees are affected by heavy competition (Fig. 12).
The simulation examples support the findings of Pukkala et al. (1994) and Vettenranta (1996), which indicate that in mixed coniferous forests it is profitable to reduce the proportion of pines during the rotation.In the present study in which the data were collected from throughout southern Finland, the growth of spruces compared with the growth of pines is more vigorous than in the previous studies of conifer mixtures in which the data were collected only from North Karelia.Thus, with the present models the optimal proportion of spruce in the end of rotation is greater than that in previous studies in which the total volume yield was used as a decision criterion.However, in all these studies the profitability of mixture is evident, as in the studies of Jonsson (1962) and Kerr et al. (1992).
The climatic variation in diameter and height growth was not taken into account, but the wide temporal range of the study material was supposed to eliminate the bias caused by the variation in growth between different years.Another deficiency of the present simulation models is the structure of the height-growth models, which are not affected by the crown ratio.The vitality of a tree, described with the aid of the crown ratio, may affect the height growth of the tree as well as it affects the diameter growth.On the other hand, for the sake of simplicity the number of models was restricted.
To find the optimal treatment schedules for each stand with a conifer mixture, the age and the structure of the stand must be determined correctly.It is possible that there is no presupposition to express common treatment instructions for conifer mixtures, but the present models are convenient tools when used together with the optimization algorithm to determine the optimal treatment schedule for the particular stand and situation.

Fig. 1 .
Fig. 1.Location of the study plots and distance from the equator.

Fig. 2 .
Fig. 2. Distribution of tree species in the study plots used in the study material.Proportion of spruces refers the relative stand basal area of spruces in the INKA-data and the relative number of spruces in the TINKA-data.

Fig. 3 .
Fig. 3.The principle of computing the competition indices.Competition index is sum of exceeding parts (E) of trees.D ij is distance between trees (m).bh i is the starting level of the paraboloid surface (see Eq. 2).S is subject tree, C is competitor and NC is non-competitor.

Fig. 5 .
Fig. 5. Effect of the stand basal area and the relative number of surrounding spruces on the crown ratio of pines and spruces.The spatial distribution is regular and all trees were presumed to have equal dimensions, h = 13 m, d = 13 cm, Nc = 6900 km and site type is fresh mineral soil.

Fig. 6 .
Fig. 6.Mean diameter and mean height of pine (solid line) and spruce (dashed line) in pure pine (Myrtillus site) and spruce (Oxalis-Myrtillus site) stands according to Koivisto (1959) and in the mixed stand (Myrtillus site) in accordance with the present study material (fine lines, respectively).Heavy line refers to stand dominant height.Stand was thinned from below and diameter growth was calculated by growth models with crown ratio.

Fig. 7 .
Fig. 7. Difference in the mean diameters received withthe different simulation models (Dg without crown ratio -Dg with crown ratio) when the stand was thinned from below (fine lines) and from above (heavy lines), so that all trees were removed from the largest diameter class.Solid lines refer to pine and dashed lines refer to spruce.

Fig. 8 .
Fig. 8.Total yield in pure pine (Myrtillus site) and spruce (Oxalis-Myrtillus site) stands (fine lines) according to Koivisto (1959) and in the mixed stand (Myrtillus site) in accordance with the present study material (heavy line).Stands were thinned from below and diameter growth was calculated by the growth models with the crown ratio.Solid line refers to pine and dashed line refers to spruce.

Fig. 9 .
Fig. 9.Total yield (m 3 ha -1 (90a) -1 ) of the stand with different species composition, thinning method and each simulation model.Thinnings are conducted at the ages of 50 and 70 years with corresponding the pre-thinning basal areas at about 30 m 2 ha -1 and 35 m 2 ha -1 and the post-thinning basal areas about 20 m 2 ha -1 and 22 m 2 ha -1 .Spruce-% refers to the relative stand basal area of spruces at the age of 90 years.

Fig. 11 .
Fig. 11.Development of the density and stand basal area during rotation without treatments.The heavy line indicates the total density of the stand, solid line indicates to the proportion of pines and dashed line the proportion of spruces.The heavy dashed line denotes the total stand basal area.

Fig. 12 .
Fig. 12.Effect of competition (C s 9) and tree height on the height growth of spruce.Tg is 15 years and height of the tree is 2, 5 and 8 meters.

Table 1 .
Xm j , X(m + 1) jk ,..., Xn jk , C(p) ijk , C(s) ijk , d ijk , h ijk , cr ijk ) + p j + t ij + e ijk(1) Mean, standard deviation, range and number of observations (N) of some variables used in modelling.d = diameter, h = height, cr = crown ratio, i d = future 5-year diameter growth, i h = future 5-year height growth, T g = age of a stand at the breast height measured from the dominant trees, G = basal area of a stand, H d = dominant height of a stand, Ec = distance from certain meridian and Nc = distance from the equator.Numbers in the brackets denote the values related to a tree instead of a stand.

Table 2 .
Competition indices used in models, parameters a and b in formula 1 and mean and maximum values of indices in each data set.Subscript x refers to the species of competitors.

Table 4 .
Residual mean square errors divided into between-stand, between-tree and within-tree components, relative RMSE and degree of determination in fixed part of each model.

Table 2 ,
H d is the dominant height of a stand and h j is the height of a subject tree.

Table 5 .
Stand variables in the stand used in the simulations, where SP = Scots pine, NS = Norway spruce, G = stand basal area, N = number of stems, V = volume, Tg = stand age at breast height, Hg = basal area weighted mean height, Hd = dominant height of the stand and Dg = basal area weighted mean diameter.Standard deviation is in the brackets.