Seemingly Unrelated Mixed-E ﬀ ects Biomass Models for Young Silver Birch Stands on Post-Agricultural Lands

: Secondary succession that occurs on abandoned farmlands is an important source of biomass carbon stocks. Both direct and indirect tree biomass estimation methods are applied on forest lands. Using empirical data from 148 uprooted trees, we developed a seemingly unrelated mixed-e ﬀ ects models system for the young silver birch that grows on post agricultural lands in central Poland. Tree height, biomass of stem, branches, foliage, and roots are used as dependent variables; the diameter at breast height is used as the independent variable. During model elaboration we used restricted cubic spline: 5 knots at the quantiles (0.05, 0.275, 0.5, 0.725, and 0.95) of diameter at breast height provided su ﬃ ciently ﬂexible curves for all biomass components. In this study, we demonstrate the use of the model system through cross-model calibration of the biomass component model using tree height measured from 0, 2, 3, and 4 available extreme trees feature in the plot in question. A di ﬀ erent number of extreme trees were measured for ﬁnal model system and our results indicated that for all analyzed components, random-e ﬀ ect predictions are characterized by higher accuracy than ﬁxed-e ﬀ ects predictions. Abstract: Secondary succession that occurs on abandoned farmlands is an important source of biomass carbon stocks. Both direct and indirect tree biomass estimation methods are applied on forest lands. Using empirical data from 148 uprooted trees, we developed a seemingly unrelated mixed-effects models system for the young silver birch that grows on post agricultural lands in central Poland. Tree height, biomass of stem, branches, foliage, and roots are used as dependent variables; the diameter at breast height is used as the independent variable. During model elaboration we used restricted cubic spline: 5 knots at the quantiles (0.05, 0.275, 0.5, 0.725, and 0.95) of diameter at breast height provided sufficiently flexible curves for all biomass components. In this study, we demonstrate the use of the model system through cross-model calibration of the biomass component model using tree height measured from 0, 2, 3, and 4 available extreme trees feature in the plot in question. A different number of extreme trees were measured for final model system and our results indicated that for all analyzed components, random-effect predictions are characterized by higher accuracy than fixed-effects predictions.


Introduction
Issues related to climate change have become increasingly important, as evidenced by the world's strongest climate-energy policy agreed upon by the European Council in October 2014. The policy's objective assumes that within the European Union (EU): (i) at least 32% of the energy demand is being covered by renewable energy sources and (ii) greenhouse gas emissions (from 1990 levels) will be cut by 40% by 2030. In this context, due to the area occupied and the amount of accumulated carbon, forest ecosystems play an important role. For example, the content of carbon in biomass of Polish forests has been estimated at 822 million tons [1]. It is worth noting, however, that due to the difficulty of defining the current method of land management, almost 800 thousand ha of forests are not included in the official statistics [2]. One of the types of land that is not covered by official statistics is secondary succession, which occurs on abandoned farmlands and are an important source of biomass carbon stocks [3][4][5]. These areas are frequently subjected to appearance of pioneer forest tree species, such as silver birch (Betula pendula Roth.) [6,7]. Those stands during the first eight years produce 31.2 Mg/ha aboveground biomass [7] and after 15 years can produce up to 75 tons Mg/ha [8]. By introducing rapid rotation in these areas, it is possible to increase the potential of fast growing silver birch stands, which can be assessed using a life cycle analysis [9]. additional transformation of the variable. The function's flexibility is controlled by the number of knots, which can be placed at certain quantiles of the independent variable [35].
One benefit from the mixed-effects model is the flexibility in prediction. If no local information to predict the random effects for the plot in question are available, a prediction using the fixed part of the model can be used. Using the fixed part of mixed-effects models provides predictions for a typical plot and in case of nonlinear models it could be biased [42]. However, if even a small amount of local measurements are available, a more accurate random-effect prediction can be used [43,44]. In the multivariate modeling context, such prediction utilizes the cross-model correlations of the model system, thus taking into account the logical interactions between the biomass components. Furthermore, those model systems allow use of the cross-model residuals correlation to predict residuals for all dependent variables for the calibration trees in question. Lappi [43] developed a bivariate system of seemingly unrelated mixed-effects model for tree height and volume, and used the system for prediction of standing tree volume both when measurements of tree height were not available and when they were available from the plot in question. In this study, we extend this idea to the multivariate system including models for biomass components and tree height. The model system allows prediction of plot-level biomass-diameter models and further improving that prediction using local measurements of tree heights from the plot in question. Above predictions is based on the estimated best linear unbiased predictor (EBLUP) and it utilizes the estimated variance-covariance structure of random effects and residuals in a way that leads to statistically optimal prediction [42].
Using empirical data from 148 uprooted trees, we developed a seemingly unrelated mixed-effects models system for the young silver birch stands growing on post agricultural lands in central Poland. Tree height, biomass of stem, branches, foliage, and roots are used as dependent variables while diameter at breast height are used as the independent variable. We demonstrate the use of the model system through cross-model calibration of the biomass component models using tree height.

Study Sites
We used empirical data from 16 temporary sample plots located in silver birch stands growing on post-agricultural lands in Mazowieckie province, central Poland (52 • 21 -51 • 24 N, 20 • 39 -21 • 26 E, Figure 1). We analyzed stands originated from natural regeneration that started after farming use was stopped. No silvicultural treatments had been performed by the time the sample material was collected. All sample plots were located in a zone of transition from the maritime to continental type within the temperate climate [45]. Mean annual temperatures were in the range 6-8 • C. The average annual precipitation in the region amounted to 550-600 mm. Soils are generally poor in nutrients. Elevation ranged from 95 to 160 m above sea level. Each sample plot consisted of about 200 trees. We measured diameter at breast height of all living trees on sample plot. At each sample plot, we randomly chose ten trees from the range of diameters. For this study we selected trees taller than 1.3 m. Age of those trees varied from 3 to 16 years (Table 1). Table 1. Characteristics of sample plots: Age (years), plot area (Area (ha)), number of trees per hectare (N), Basal area per hectare (BA (m 2 ha −1 )), and trees: diameter at breast height (DBH (cm)), height (H (m)), stem (ST), branches (BR), foliage (FL), and roots (RT) dry biomass (kg) of the sampled trees.

Uprooted Trees
The sample trees were dug out with their root systems. We marked the base level at the bark prior to excavation to separate aboveground part of the tree from belowground one. We divided each tree into the following biomass components: stem, branches, foliage and roots (>2 mm of diameter according to [46]) separated from the stem at the previously marked base level.
We cleared the roots from soil and other organic materials using brushed and compressed air. We weighted all biomass components in the field using precise portable scales (accuracy 0.01 g). We took samples of each component (stem discs in the middle of each 1 m section, random samples of branches, foliage, and roots) from every tree to calculate the ratio between fresh and dry biomass. We dried the samples at 105 °C [47] until they reached a constant weight. We calculated the dry biomass of each biomass component for each tree on the basis of corresponding fresh to dry biomass ratios [46]. Finally, empirical material consisted of the data from 148 trees taller than 1.3 m. Diameter at breast height of those trees varied from 0.1 to 9.7 cm (Table 1, Figure 2).

Uprooted Trees
The sample trees were dug out with their root systems. We marked the base level at the bark prior to excavation to separate aboveground part of the tree from belowground one. We divided each tree into the following biomass components: stem, branches, foliage and roots (>2 mm of diameter according to [46]) separated from the stem at the previously marked base level.
We cleared the roots from soil and other organic materials using brushed and compressed air. We weighted all biomass components in the field using precise portable scales (accuracy 0.01 g). We took samples of each component (stem discs in the middle of each 1 m section, random samples of branches, foliage, and roots) from every tree to calculate the ratio between fresh and dry biomass. We dried the samples at 105 • C [47] until they reached a constant weight. We calculated the dry biomass of each biomass component for each tree on the basis of corresponding fresh to dry biomass ratios [46]. Finally, empirical material consisted of the data from 148 trees taller than 1.3 m. Diameter at breast height of those trees varied from 0.1 to 9.7 cm (Table 1, Figure 2).

Uprooted Trees
The sample trees were dug out with their root systems. We marked the base level at the bark prior to excavation to separate aboveground part of the tree from belowground one. We divided each tree into the following biomass components: stem, branches, foliage and roots (>2 mm of diameter according to [46]) separated from the stem at the previously marked base level.
We cleared the roots from soil and other organic materials using brushed and compressed air. We weighted all biomass components in the field using precise portable scales (accuracy 0.01 g). We took samples of each component (stem discs in the middle of each 1 m section, random samples of branches, foliage, and roots) from every tree to calculate the ratio between fresh and dry biomass. We dried the samples at 105 °C [47] until they reached a constant weight. We calculated the dry biomass of each biomass component for each tree on the basis of corresponding fresh to dry biomass ratios [46]. Finally, empirical material consisted of the data from 148 trees taller than 1.3 m. Diameter at breast height of those trees varied from 0.1 to 9.7 cm (Table 1, Figure 2).

Separate Mixed-Effect Models
We started the analysis by modeling the relationship between the dependent variables and tree diameter separately for each biomass component using a liner mixed-effect model. In our data, sample plots defined a single grouping factor and no other potential factors of grouping were present. Each model for the tree j, j = 1,...,ni, of plot i, i = 1,...,M, was of form:

Separate Mixed-Effect Models
We started the analysis by modeling the relationship between the dependent variables and tree diameter separately for each biomass component using a liner mixed-effect model. In our data, sample plots defined a single grouping factor and no other potential factors of grouping were present. Each model for the tree j, j = 1, . . . , n i , of plot i, i = 1, . . . , M, was of form: where x ij β is the fixed part and z ij b i + ij the random part. We assume that: where b i are independent among plots and of the residual errors ij . Matrix ψ is a positive definite variance-covariance matrix of random effects. The zero-mean residual errors ij were mutually independent, normal with variance.
where σ 2 and δ are the estimated scale and shape parameters for the power function, respectively. The special case δ = 0 leads to homoscedastic residuals. Notice that our model includes only diameter at breast height (DBH) as a fixed predictor, therefore prediction does not require that e.g. measurement of tree height is available. However, our multivariate system allows also improved predictions by utilizing measurements of height as well.
To ensure that the predicted biomass components and tree heights are positive, logarithmic biomass and logarithm of height-1.3 was used as y ij . However, the logarithmic biomasses were not linear with respect of DBH or its logarithm and we used a restricted cubic spline to model this curvilinearity. The following restricted cubic spline with K nodes t k , (min(DBH) < t 1 < . . . < t K < max(DBH)) was used in the fixed part [41].
The predictors in the last K-2 terms are the truncated power basis functions of form: where (DBH − t k ) 3 + is zero when DBH < t k and (DBH − t k ) 3 otherwise. This solution allows us to define model, which is smooth and curvilinear with respect to DBH, but linear with respect to the regression coefficients; equation (5) provides a way to compute nonlinear transformations of DBH. Initial exploration of the data suggested random intercept and coefficient of diameter in all models, therefore the random effects of form: i DBH ij were used in all models. We used five knots placed at certain quantiles of the DBH following the suggestion of Harrell [41]. The number of knots was chosen based on the residuals behavior analysis and AIC (Akaike Information Criterion). Residual analysis for all dependent variables was carried out using plots of Pearson residuals in three variants: (i) the means of residuals in 10 classes of the standardized DBH were plotted together with the 95% confidence interval of individual observations (mean ± 1.96 SD) and the confidence interval for class mean (mean ± 1.96 SE). Standardized diameter was calculated as difference between individual tree diameter and mean diameter of the sample plot, divided by standard deviation of diameter in sample plot. The similar graphs were produced also (ii) in relation to DBH and (iii) the fitted values. Residuals were plotted as implemented in function mywhiskers of R-package lmfor [48].

Multivariate Seemingly Unrelated Mixed-Effects Model System
Multivariate (seemingly unrelated) mixed-effects models are seemingly unrelated models with random group effects. They are useful, in particular, when several independent variables are measured from the same sample units and the relationship among the random effects and residuals across models is of interest. Assume that we have L individual mixed-effects models and pool the observations of the i . The multivariate seemingly unrelated mixed-effects model system for each plot i, i = 1, . . . , m as [35]: y We assume that the cross-model correlations of the random effects and residuals are the same for all groups. By defining and correspondingly for b i and i (see [35] for details), the multivariate model can be written as univariate model: where X i and Z i are block-diagonal matrices, which have the response-specific model matrices on the diagonal. The variance-covariance matrix of the random effects is and the correlation matrix of residual errors is are the cross-model covariances of random effects and ρ lm is the correlation between the residual errors for the same tree between responses l and m. The variance-covariance matrix of residual errors is where the diagonal matrix Γ includes blocks var The definition of the model for whole data is straightforward. The parameter estimation was done using REML using R-package nlme. See [35] for more details about the model formulation and for an example about the implementation in R.

Fixed and Cross-Model Random-Effect Prediction
In general, mixed-effect models allows to use fixed-effect and random-effect prediction. Fixed-effect prediction assumes that the expected value of the random effects is 0 and it is optimal when no sample information about dependent variables are available from the plot in question. If measurements of even one tree are available from the plot in question, the random-effect prediction is justified. It may sound surprising that one measurements can be used for prediction of several parameters, but this is possible and justified because the approach uses the estimated variances within plots and between plots to optimally weight the fixed part of the model and the local measurements. The multivariate seemingly unrelated mixed-effects model allows to use also the cross-model correlations of random effects and residual errors to predict all dependent variables using measurements of some of them. Generally speaking, we can predict the random group effects for all dependent variables: both for the variables that have local measurements and for the variables that do not have local measurements for the plot in question [35]. We could also predict the residual errors of all response variables for the sampled trees [35,42].
Technically, we define vector y 0 that includes from y i only those elements that have been observed. We define matrices D 0 and R 0 by removing from D and R those rows and columns that correspond to the unobserved variables. In addition, we define matrix C by removing from D all the columns that correspond to the unobserved variables but keeping all the rows, and Z 0 and X 0 by removing from Z i and X i those rows and columns that correspond to the unobserved variables. The best linear predictor of the complete random effect vector for group i becomes [35]: The prediction variance is The fixed part and variance-covariance matrices are replaced by their estimates, which leads to EBLUP. In our case, the heights and biomass components were predicted in a logarithmic scale. To transform them to the original scale we applied bias correction by adding half of the prediction variance into the predictions before applying the exponential transformation [35,49]. The required prediction variance was taken from the diagonal of Z i var b i − b i Z i + var(e i ). That variance ignores the components related to the estimation errors of fixed effects and variance-covariance components, but their effect is marginal because of the rather large data set.
In our case during plot-specific cross-model random-effect prediction we used tree height as the observed dependent variable that is used for calibration. We used a different number of available measured trees for prediction. The fixed-effect approach was applied for 0 available trees, whereas the cross-model random-effect prediction was created separately for 2, 3, and 4 available trees. From the literature of experimental design [50], we know the prediction error variance of a linear regression model is the lowest if the x-variables are taken from the endpoints of the range (extreme trees). Therefore, we used this procedure for cross-model random-effect prediction. We first randomly divided each sample plot, with 10 sample trees, into prediction and evaluation part. Afterwards, based on prediction part we conducted extreme trees sampling (2, 3, 4) for height measurement. We used DBH as the selection criterion as follows: 2 trees-the thinnest and the thickest tree; 3 trees-1 the thinnest and 2 the thickest trees; 4 trees-2 the thinnest and 2 the thickest trees. The sampling procedure was replicated 20 times for each plot. Based on trees forming the evaluation part of each sample plot and which were not used in calibration, we computed the root mean square error (RMSE). RMSE was calculated for each dependent variable as mean error achieved during all replication.

Results
Analysis of the residuals for individual models for each depended variable indicated the presence of heteroscedasticity. Therefore, we modeled the residual variance using a power function (3). The assumed variance function modeled the heteroscedasticity well ( Figure 3). Likewise, a likelihood testing indicated that the variance function significantly improved the fit (p < 0.0001). In the spline regression (4), we found that using five knots at the quantiles (0.05, 0.275, 0.5, 0.725, 0.95) of DBH provided sufficiently flexible curves for tree height and all analyzed biomass components ( Figure 4).
Once the model form for the individual models was found, the model was fitted as a five-component seemingly unrelated mixed-effects model. The estimated regression coefficients and variance function parameters of the model system are given in Table 2.  The cross-model correlations of the random effects range from 0.08 to 0.995, being moderate or strong in most cases, which indicates that the cross-model calibration will be successful ( Table 3). The random effects for height are most strongly correlated with random effects estimated for stem biomass. Considering the analyzed biomass components, it should be noted that the strongest relationship is observed between random effects achieved for foliage and branches biomass. Random effects for stem biomass, as component with the largest share in aboveground biomass, are strongly correlated with foliage and roots biomass. While random effects for roots biomass is strongly correlated with branches biomass.
Correlation between residual errors for analyzed dependent variables varies depending on the feature being assessed ( Table 4). The strongest correlation is observed in the case of branches and foliage biomass. In addition, residual errors for roots biomass has a strong relationship to other aboveground biomass components. Residual error for tree height are not strongly correlated with biomass components-the highest in case of stem biomass. Table 3. Random effects variance covariance matrix for model system. Matrix diagonal-variance, matrix upper triangle-correlation between random effects.

Fixed and Cross-Model Random Effects Prediction
Applied different number of extreme trees measured for final model system indicated that for all analyzed depended variables random effects predictions have lower RMSE then fixed effects prediction, which serves as a reference in Figure 5 (Number of trees = 0). When the number of available trees is increased, the error decreases for all depended variables. The differences between RMSE values showed in Figure 5 would be clearer with more trees available. However, due to the limitations in the number of trees per plot of the data, we could not include more trees for random effect prediction.
In the case of height, using four trees in the random-effects prediction decreased the RMSE by 0.045 m (5.4%, Table 5) compared to the fixed-effect prediction. In case of biomass components the highest relative decrease was reached for foliage biomass and equals 41.1%, whereas the lowest relative difference occurs in case of stem biomass (6.7%; Table 5).

Discussion and Conclusions
In this study, based on 148 uprooted trees from 16 temporary sample plots, we developed a seemingly unrelated mixed-effects biomass model system for young silver birch stands growing on post agricultural lands in central Poland. We defined diameter at breast height as independent variable and tree height as well as the whole tree biomass components such as stem, branches, foliage, and roots biomass as dependent variables. During our analysis we applied fixed and random effect prediction. Under prediction we included the sample plot as the grouping level and tested height measurements for 2, 3, and 4 available trees in plot in question. It should be noted that the resulting model system is a solution with high potential in case of biomass modeling. On the one hand, it allows us to achieve information about roots biomass, which is difficult to measure, but is crucial according to the understanding and estimating carbon accumulation in birch forest ecosystems [8,23,24]. On the other hand, it potentially allows prediction of the tree-level residuals and therefore the use of either diameter or diameter and height in prediction of tree biomass using the procedures explained in Lappi [43] (see also [42]). Moreover, it allows us to take into account various data sources [55] and is a new voice in the discussion on biomass prediction using diameter at breast height and height [13,15]. In the context of potential use of achieved solution it is worth noting that during the cross-model random-effect prediction we focused on the trees available only for height measurements as the most accessible dependent variable. Therefore, taking into account the larger number of available trees and other dependent variables (such as measurements of stem biomass to improve the prediction of other biomass components) achieved model system will allow a broader understanding of accuracy of cross-model random-effect prediction.
In this study, we estimated total biomass as the sum of the individually calculated best regression models of the biomass of its components [29]. The variance is estimated by taking into account the cross-model covariances of random effects and residual errors. This approach was also used during creation of Canadian national aboveground biomass equations for six tree species [56]. However, they defined two types of models: (i) the DBH-based equations where, as in our case, DBH is an independent variable and (ii) the DBH-and height-based equations where, unlike in our study, the tree height was included as independent variable. Aboveground biomass additive models based on SUR were also created for the eight relevant tree species in Germany [57]. Likewise, Carvalho and Parresol [58] applied nonlinear seemingly unrelated regression (as nonlinear joint-generalized regression) for stem and crown biomass functions for Pyrenean oak trees in Portugal. Morever, in Poland, research was carried out taking into account SUR in biomass modeling. Zasada et al. [34] and Bronisz and Zasada [15] created additive biomass models for aboveground biomass of Scots Pine. While Bronisz et al. [22] focused on empirical equations for estimating aboveground biomass of silver birch. The aforementioned studies are examples of creating comprehensive solutions enabling biomass prediction of various tree components. In addition, these are solutions where DBH and height are defined as independent variables, while in our work we used a different approach which, to our knowledge, has not been previously used in biomass modeling. Parresol [59] formulated a system of biomass models where one of the components is the total biomass, which is just the sum of the other components of the system. Compatibility of predictions was ensured through parameter restrictions. However, even though the estimated variance of the total biomass in the empirical data set was lower, no theoretical arguments for the better performance of such an approach over the one we used was provided.
It is worth adding that model systems are not only used during biomass modeling. Fu et al. [60] indicates that seemingly unrelated regression allows to reach higher prediction accuracy for crown width than other methods such as adjustment in proportion or ordinary least square regression. Mehtätalo [61] created a model system for stand-specific diameter distribution while Siipilehto [62] developed a system of models for 20 stand-level variables (e.g., basal area, mean diameter and dominant height). Instead Maltamo et al. [55] predicted and calibrated diameter at breast height, tree height, volume, dead branch height and crown base height using airborne laser scanning.
Regardless of whether biomass is modeled using individual formulas or a model system, solutions based on linear regression are possible [63]. However, due to modeled biomass relationships, nonlinear regression is more often considered [21,22,25]. One of the reasons is that the biomass is always positive, which was in our case ensured by modeling logarithmic biomass. Another acceptable solution would be a generalized linear mixed-effect model based on gamma distribution [64]. During our study we modeled the curvilinear relationship of biomass and diameter by using a spline regression, which ensures that the model remains as linear. It is continuous function consisting of pieces of second or third order polynomials that are joined at knots [41]. In addition, Harrell [41] indicates that number of knots depends on the sample size. Others claim that for many datasets four knots offers a good enough fit of the model and is a compromise between flexibility and loss of precision caused by overfitting a small sample. Our results for individual mixed-effect models indicated that in the case of biomass modeling of young birch stable residual behavior and lowest AIC value is obtained for five knots. The results obtained in our study indicate that it is possible to use the spline function for biomass modeling based on empirical data. Splines have been used in biomass modeling previously in Bollandsås et al. [65]. During this research authors tested three approached for biomass estimation based on field and laser data collected in a mountain forest in southeastern Norway.
The data for biomass models is often characterized by a grouped structure. It usually contains information about individual tree biomass components for different sample plots. The relationship between variables varies among different sample plots. The resulting dependence is properly taken into account by using the mixed-effects models approach [66]. Smith et al. [25] in research on belowground and whole tree biomass of birch in Norway indicates that nonlinear mixed effects models approach for those type of modeling, with similar as in our case, single grouping factor is useful methodology for those type of modeling. The best model for whole tree biomass was evaluated using two independent variables (diameter at breast height and height). Repola [63] created three separate multivariate variance component linear models for aboveground and belowground biomass for birch trees in Finland. To reach homoscedasticity of the variance and to ensure positive biomasses, the author used logarithmic transformation of data. In our case, the logarithmic transformation did not lead to homoscedastic errors and an additional variance function was therefore used. As opposed to our approach, another author defined diameter at breast height and height as independent variables.
Linear mixed-effect models were also tested during the Fehrmann et al. [36] study. In their research, aboveground biomass was estimated for Finnish Scots pine and Norway spruce trees using three approaches: (i) k-nearest neighbour (k-NN), (ii) plot-specific linear mixed-effect biomass models, in which DBH and height were used as independent variables, and (iii) simple linear models. Based on the achieved results, they claimed that even fixed effect prediction based on linear mixed models allows for greater precision than an ordinary least squares regression.
In addition to taking into account the dependence structure caused by the grouping of the data, the mixed-effects modeling allows for the flexible prediction alternatives. In cases where in there is a lack of measurements, only fixed-effects can be used, yet additional measurements allow random-effects to be predicted [66]. Consequently, it allows for greater predictions accuracy in plots where additional measurements are available [44]. However, in the case of biomass modeling, obtaining additional measurements is troublesome because it is associated with the felling of trees. Multivariate mixed-effects model application allows to take full advantage of the mixed-effects models. In our case, additional measurements for tree height allows to reach more accurate random cross-model prediction for all analyzed biomass components. The same idea could be generalized also to root biomass: if measurements of the above-ground biomass components are available, our model would allow for (1) more accurate prediction for all biomass components of the plot and (2) even more accurate prediction of the root biomass for the sample trees through prediction of the residual error of that tree [42,43].
Author Contributions: K.B., conceptualization, software, investigation, resources, data curation, and writing-original draft. L.M., methodology, software, validation, writing-review and editing. All authors have read and agreed to the published version of the manuscript. Abstract: Secondary succession that occurs on abandoned farmlands is an important source of biomass carbon stocks. Both direct and indirect tree biomass estimation methods are applied on forest lands. Using empirical data from 148 uprooted trees, we developed a seemingly unrelated mixed-effects models system for the young silver birch that grows on post agricultural lands in central Poland. Tree height, biomass of stem, branches, foliage, and roots are used as dependent variables; the diameter at breast height is used as the independent variable. During model elaboration we used restricted cubic spline: 5 knots at the quantiles (0.05, 0.275, 0.5, 0.725, and 0.95) of diameter at breast height provided sufficiently flexible curves for all biomass components. In this study, we demonstrate the use of the model system through cross-model calibration of the biomass component model using tree height measured from 0, 2, 3, and 4 available extreme trees feature in the plot in question. A different number of extreme trees were measured for final model system and our results indicated that for all analyzed components, random-effect predictions are characterized by higher accuracy than fixed-effects predictions.

Conflicts of Interest:
The authors declare no conflict of interest.