Aggregated Biomass Model Systems and Carbon Concentration Variations for Tree Carbon Quantiﬁcation of Natural Mongolian Oak in Northeast China

: Three systems of additive biomass models were developed and the e ﬀ ects of tree components, tree sizes, and tree growing regions on the carbon concentration were analyzed for Mongolian oak ( Quercusmongolica Fisch. exLedeb.) inthenaturalforestsofNortheasternChina. Thenonlinear seemingly unrelated regression (NSUR) method was used to ﬁt each of the three systems simultaneously; namely, aggregated model systems with no parameter restriction (AMS0), aggregated model systems with one parameter restriction (AMS1), and aggregated model systems with three parameter restrictions (AMS3). A unique weighting function for each biomass model was applied to address the heteroscedasticity issue. The systems assertively guarantee the additivity property, in which, the summation of the respective predicted tree components (i.e., root, stem, branch, and foliage) will match the prediction of subtotals (i.e., crown and aboveground) and total biomass. Using one-, two-, and three-predictor combinations (i.e., D (diameter at breast height), D and H (total height), and D, H, and CL (crown length)) as the general model underlying formats, three systems of additive biomass model were developed. Our results indicate that (1) all of the aggregated model systems performed well and the di ﬀ erences between the systems were relatively small; (2) the rank order of the three systems based on an array of statistics are as follows: AMS0 > AMS1 > AMS3; (3) the carbon concentration signiﬁcantly varied depending on the types of tree tissues and growing regions; (4) the regional respective component carbon concentration and regional weighted mean carbon concentration multiplied by observed biomass value appeared to be the best approach to calculate carbon stock.


Introduction
On a global scale, forests are known as one of the most dominant terrestrial ecosystems, covering 13 billion hectares or approximately 30.6% of the Earth's total land area [1]. Forests contribute in the storage of a large amount of biomass and play an important role in maintaining atmospheric carbon balance on Earth [2]. Researchers have been facing a challenge for decades to enhance the accuracy of tree biomass and carbon estimations, which are crucial for understanding climate change, nutrient cycling, forest productivity, forest health, etc. [3][4][5]. Measuring the actual weight of each tree component is unquestionably known as the most precise approach for estimating tree biomass. not additive. Therefore, the aims of this research were: (1) to construct three aggregated biomass model systems (i.e., AMS0, AMS2, and AMS3) with weighted NSUR based on different combinations of independent variables for Mongolian oak in northeast China; (2) to verify the predictive performance of the three systems using jackknifing techniques; (3) to evaluate six procedures of quantifying carbon stock of an individual tree; and (4) to compare the aggregated model systems developed in this study against the previously developed biomass equations.

Data Collection
This research was undertaken in natural secondary forests located in Heilongjiang Province (from 121 • 11 E to 135 • 05 E and from 43 • 25 N to 53 • 33 N) and Jilin Province (from 121 • 38 E to 131 • 19 E and from 40 • 50 N to 46 • 19 N), the two largest provinces in Northeast China (Figure 1). This area encompasses five regions investigated in this study; namely, the Northern slope of Xiaoxing'an Mountain (XXMN), the Southern slope of Xiaoxing'an Mountain (XXMS), the Eastern slope of Zhangguangcai Mountain (ZGCME), the Western slope of Zhangguangcai Mountain (ZGCMW), and the Changbai Mountain regions (CBM). These regions are characterized as having a continental monsoon climate within the Dwa, Dwb, and Dwc categories of the Köppen-Geiger climate classification system, having dry winters along with hot, warm, and cold summers, respectively [38]. The soils are mostly either Eutroboralfs or Haplumbrepts across the five regions. The elevation, mean annual precipitation, and mean annual temperature varies from 300 to 1500 m above sea level, between 500 and 800 mm, and from −4 to +6 • C, respectively.
A total of 72 trees from natural secondary forests were destructively sampled to measure the following parameters: (1) tree basic variables such as four directional crown width (CW), diameter at breast height (D), the length of the live crown (CL) and the total tree height (H); (2) the actual fresh weight of the roots, stems, branches, and foliage; and (3) the moisture content and carbon concentration of the sampled roots, stems, branches, and foliage. Note that this study did not separately measure tree bark; thus, all of the stem weights and measurements were of trees not debarked. A detailed explanation of the determination of biomass and carbon concentration, including field and laboratory measurements, can be found in Dong et al. [7,39]. The summary statistics of all sampled trees are listed in Table 1. concentration in biomass tissues for Mongolian oak in northeast China. However, both of these studies used a relatively small sample size, and particularly the biomass models developed by Wang [36] were not additive. Therefore, the aims of this research were: (1) to construct three aggregated biomass model systems (i.e., AMS0, AMS2, and AMS3) with weighted NSUR based on different combinations of independent variables for Mongolian oak in northeast China; (2) to verify the predictive performance of the three systems using jackknifing techniques; (3) to evaluate six procedures of quantifying carbon stock of an individual tree; and (4) to compare the aggregated model systems developed in this study against the previously developed biomass equations.

Data Collection
This research was undertaken in natural secondary forests located in Heilongjiang Province (from 121°11′E to 135°05′E and from 43°25′N to 53°33′N) and Jilin Province (from 121°38′E to 131°19′E and from 40°50′N to 46°19′N), the two largest provinces in Northeast China (Figure 1). This area encompasses five regions investigated in this study; namely, the Northern slope of Xiaoxing'an Mountain (XXMN), the Southern slope of Xiaoxing'an Mountain (XXMS), the Eastern slope of Zhangguangcai Mountain (ZGCME), the Western slope of Zhangguangcai Mountain (ZGCMW), and the Changbai Mountain regions (CBM). These regions are characterized as having a continental monsoon climate within the Dwa, Dwb, and Dwc categories of the Koppen-Geiger climate classification system, having dry winters along with hot, warm, and cold summers, respectively [38]. The soils are mostly either Eutroboralfs or Haplumbrepts across the five regions. The elevation, mean annual precipitation, and mean annual temperature varies from 300 to 1500 m above sea level, between 500 and 800 mm, and from -4 to +6°C, respectively.
A total of 72 trees from natural secondary forests were destructively sampled to measure the following parameters: (1) tree basic variables such as four directional crown width (CW), diameter at breast height (D), the length of the live crown (CL) and the total tree height (H); (2) the actual fresh weight of the roots, stems, branches, and foliage; and (3) the moisture content and carbon concentration of the sampled roots, stems, branches, and foliage. Note that this study did not separately measure tree bark; thus, all of the stem weights and measurements were of trees not debarked. A detailed explanation of the determination of biomass and carbon concentration, including field and laboratory measurements, can be found in Dong et al. [7,39]. The summary statistics of all sampled trees are listed in Table 1.  The allometry relationship between tree biomass components and several available predictors were primarily screened through visual inspection to ascertain whether or not each component biomass model could appropriately be developed as multivariate allometric equations ( Figure 2).   The allometry relationship between tree biomass components and several available predictors were primarily screened through visual inspection to ascertain whether or not each component biomass model could appropriately be developed as multivariate allometric equations  In addition, the regional effect on biomass data for each tree component was tested by using a dummy variable method, and the difference was found to be insignificant. Thus, the effect of different tree growing regions was not considered. Nonlinear regression with an additive error term was determined as the general biomass equation defined as follows: where W i represent each biomass components in kg (i = r for roots, s for stems, b for branches, and f for foliage); X j are several available predictors (j = 1, . . . , k), such as D, H, CW, and CL; β ij are the model coefficients; ε i is the error term. Each of the tree component biomass equations were fitted using various candidate combinations of available predictor (i.e., D, H, CW, and CL). Then, three model fitting statistics were used to determine the "best model" based on different predictors for each component biomass equation: (i) Akaike information criterion (AIC), (ii) root mean squared error (RMSE), (iii) coefficient of determination (R 2 ), and (iv) significance of model parameters. The model fitting statistics showed that the component biomass equations have different combinations of independent variable. Finally, we selected one-, two-, and three-independent variables as the underlying formats to develop three systems of additive biomass model.

Aggregated Model Systems
The component biomass data were required to be simultaneously fitted in order to elucidate the intrinsic correlations among component biomass of the same individuals. In this study, we proposed three different structural arrangements for the aggregated model systems, i.e., those with no parameter restriction (AMS0), those with one parameter restriction (AMS1), and those with three parameter restrictions (AMS3), as specified below: • Aggregated model systems with no parameter restriction Referring to the model structure detailed in Affleck and Diéguez-Aranda [25], the structural system of AMS0 are as follows: where W r represents roots, W s refers to stems, W b refers to branches, and W f refers to foliage biomass. The other abbreviations were the same as in Equation (1).

•
Aggregated model systems with one parameter restriction Referring to the model structure detailed in Parresol [17], the AMS1 ensures the additivity between tree component, subtotal, and total biomass with one parameter restriction, i.e., tree total biomass equals the summation of the whole tree component biomass. The structural system of AMS1 are given as follows: where W t represents total biomass. The other abbreviations were the same as in Equations (1) and (2).

•
Aggregated model systems with three parameter restrictions Referring to another of Parresol's [17] model structure, which has been applied in Dong et al. [19], the AMS3 ensures the additivity between the tree component and total biomass with three parameter restrictions, i.e., (1) tree crown biomass is equivalent to the summation of the branch and foliage biomass; (2) tree aboveground biomass is equivalent to the summation of the stem, branch, and foliage biomass; (3) tree total biomass is equivalent to the sum of the whole tree components biomass. The structural system of AMS3 are given as follows: where W c and W a represent the crown and aboveground, respectively. The other abbreviations were the same as in Equations (1)-(3).

Weighting Function for Heteroscedasticity
Tree biomass data frequently exhibit heteroscedasticity in the model's error variance. Hence, in the context of nonlinear regression, a weighting function must be specified for each biomass model in order to counter the problem. It is well known that the model residuals of the ith observation have a functional correlation with one or several combinations of predictors. Thus, a power function of σ 2 i = (x i ) p was determined for the independent variables, where the power coefficient p can be obtained by iteratively using the error variance model of e 2 i = (x i ) p , in which, e i is the unweighted (OLS) model residual. Two different formats of weight functions (i.e., x = D and x = D·H) were compared, but the difference was not significant between them. Hence, p was determined for each biomass equation and 1/D p was chosen to be the weighting function. The programming code of resid. W i = resid.W i / √ D p was specified by utilizing the PROC MODEL procedure [40] in order to homogenizes the residuals and improve the model fitting, as noted by Balboa-Murias et al. [41] and Dong et al. [19].

Model Assessment and Validation
As described above, the three candidate aggregated model systems based on different predictor variables were evaluated in this study. The biomass dataset (sample size n) was used in its entirely to establish the biomass equations. Meanwhile, the jackknifing technique was utilized to validate the model, in which the sample size of all-except-one observation (n -1) was used to build the biomass models. Further, the fitted model was utilized to estimate the biomass value for the held-out observation. Five validation statistics of jackknifing were employed to calculate and compare the model performance: coefficient of determination (R 2 ), root mean square error (RMSE), mean prediction error (MPE), mean absolute error (MAE) and mean absolute percent error (MAPE). These statistical analyses are summarized as Forests 2020, 11, 397 where W ik , W ik , and W are the observed, predicted, and average biomass, respectively, of the ith sample tree from the k components (k = r for roots, s for stems, b for branches, and f for foliage), which were fitted by using all of the n observations; W ik,−i is the predicted value of the model which was fitted by using n -1 observations.

Effects of Tree Size, Region, and Component
To investigate the effect of tree sizes (D), growing regions, and tree components on the variation of carbon concentration, multi-factor analysis of variance (ANOVA) through a general linear model was utilized. The model form can be written as follows: where y ijk represents the value of observed carbon concentration, µis the average value, α i is the effect of the ith tree D on y ijk , β j is the effect of the jth region on y ijk , γ k is the effect of the kth tree component on y ijk , (α · β) ij and (α · γ) ik represent the interaction effect between tree D with growing region and tree component, respectively, (β · γ) jk is the interaction effect between region and components, (α · β · γ) ijk presents the interaction effect among the three variables, and ε ijk is the expression of error. The statistical analysis system (SAS) software [40] was used for analyses, utilizing the PROC GLM procedure to investigate the variation in carbon concentration.

Tree Carbon Quantifications
Six alternative procedures to quantify respective components, subtotals, and total carbon stock values were evaluated. Since each individual tree component has a different carbon concentration, it is undoubtedly very complex to always define the value of carbon stock by only relying on it. Once the carbon stock quantification for each tree component was obtained, they were aggregated to produce the subtotal and total carbon stock for a particular individual tree.

1.
Regional respective component carbon concentration (RRCCc): This procedure multiplied the observed biomass value (W i ) of each component by the respective component carbon concentration for each region. Further, the carbon stock of an individual tree was obtained by summing the component estimates.

2.
Regional weighted mean carbon concentration (RWMCc): In this procedure, each of the observed component biomass values (W i ) were multiplied by the weighted mean carbon concentration for each region. The weighted mean carbon concentration was measured as follows: where WC is the weighted mean carbon concentration for the jth region, and C ij and B ij are the measured carbon concentration and the relative proportion to total biomass partitioning for the ith component within the jth region, respectively.

3.
Respective  The weighted calculation was similar with the RWMCc procedure but disregarded region. Hence, the carbon stock of the subtotal and total individual tree was acquired by the summation of each component's carbon stock.

5.
Generic carbon concentration conversion factor 1 (GCCC f -1): This procedure multiplied the observed biomass value (W i ) of each component by a generic carbon concentration conversion factor of 0.45.

6.
Generic carbon concentration conversion factor 2 (GCCC f -2): This procedure resembles GCCC f -1 instead employing the generic carbon concentration conversion factor of 0.50 for tree components to quantify the carbon stock of an individual tree.
In order to identify the best and most reliable carbon stock quantification method, we conducted a comparison between the six procedures. The absolute relative error (RE) measured respective biases of carbon stock quantification procedures, computed as where Oc i and Qc i are the ith observed value of carbon stock calculated by multiplying the specific component biomass with its respective carbon concentration and the value of carbon stock quantified by six different procedures for the ith sampled tree, respectively; N represents the total sample size of carbon stock quantification data (i.e., N = 64).

Variation of Carbon Concentration
The average, minimum, maximum, and standard deviation of the carbon concentration values for each sampled tree components and RWMC C values from five different regions are given in Table 2. The carbon concentrations found in Mongolian oak were highly varied; values per tree ranged from 0.3949 to 0.5207 across whole components and regions. The highest carbon concentration of tree components for all regions was found in foliage (varying from 0.4293 to 0.5207), whereas the lowest carbon concentration was varied between roots (regions XXMN and XXMS) and branch (regions ZGCME, ZGCMW and CBM). In total, carbon concentration in each component can be ranked as follows: foliage > stem > branch > root. XXMS had the highest RWMC C with the amount of 0.4747 ± 0.0102 (mean ± standard deviation), while ZGCME had the lowest RWMC C with the amount of 0.4349 ± 0.0095. The types of tree components and growing regions were proven to significantly affect the carbon concentration at the 5% threshold, while tree sizes (D) had no such significant effect (Table 3). This result indicates, at least in our sampled data, that tree diameter at breast height does not significantly affect tree carbon concentration. Moreover, our findings were strengthened by the multiple comparison tests (employing Tukey-Kramer adjustment method), demonstrated by diffograms ( Figure 3).

Aggregated Model Systems and Validation of Models
The estimated parameters and jackknifing validations for the AMS0, AMS1, and AMS3 based on D only, combination of D and H, and combination of D, H, and CL are shown in Tables 4, 5, and 6, respectively. CW was found to give insignificant improvement as a predictor in all tree component equations; thus, it was excluded from our models. Preliminary inspection showed that adding H into

Aggregated Model Systems and Validation of Models
The estimated parameters and jackknifing validations for the AMS0, AMS1, and AMS3 based on D only, combination of D and H, and combination of D, H, and CL are shown in Tables 4-6, respectively. CW was found to give insignificant improvement as a predictor in all tree component equations; thus, it was excluded from our models. Preliminary inspection showed that adding H into the root and branch equations, along with the insertion of CL into the stem equation would not significantly improve the model fitting. Hence, D was used as the sole predictor for root and branch biomass models (Tables 4-6), the combination of D and H were used as predictors for the stem (Tables 5 and 6) and foliage biomass models (Table 5), and the complete set of predictors (D, H, and CL) were used in the foliage biomass models (Table 6). D was unquestionably a major predictor of biomass for the three aggregative model systems. The regular D model itself accounted for ≥93%, ≥95%, ≥96%, and ≥98% of variation in foliage, root, stem, and all (branch, crown, aboveground, and total) biomass, respectively ( Table 4). The insertion of H into both of the stem and foliage biomass models enhanced the model fitting performance quality of the stem, aboveground, and total models as a consequence of smaller RMSE, MPE, MAE, and MAPE (Table 5, Figure 4A,B), while the inclusion of CL into the foliage biomass models only gave a slight improvement for the model fitting of the foliage and crown biomass models (Table 6 and Figure 4C).
The parameter estimations of D (β i1 ) were positive in all biomass component equations, whereas the parameter estimations of H (β i2 ) were found to be positive in all except for the foliage component (Tables 5 and 6). The positive value of parameter estimations of D and H in stem biomass models implicated their positive relationship, in which, for the same amount of D, tree stem biomass increased with increasing tree height. The positive parameter estimations of D but negative parameter estimations of H in the foliage component indicated that an increase in D would increase foliage biomass. In contrast, increasing tree height would decrease foliage biomass for the same value of D.
In the present study, three systems of aggregative biomass models have been developed for Mongolian oak in natural forests in Northeast China. Each biomass component equation was jointly fitted by weighted NSUR, and the weighting functions were D 3.4640 , D 2.6952 , D 3.6937 , D 2.2455 , D 4.7602 , D 3.1438 , and D 2.8621 for root, stem, branch, foliage, crown, aboveground, and total biomass, respectively. Slight differences occurred among the parameter estimation values of AMS0, AMS1, and AMS3, and the majority of the SEs for AMS3 were lower than either AMS1 or AMS0 (Tables 4-6). The validation statistics suggested that AMS0, AMS1 and AMS3 had very close values of R 2 , RMSE, MAE and MAPE for predicting component, subtotal, and total biomass. However, the AMS0 yielded smaller MPE than the other two aggregative biomass modeling approaches. Hence, the MPE, MAE, and MAPE across diameter classes were visualized in order to enhance the quality of validating tree biomass models (Figure 4). Results show that the three aggregative systems based on D, D and H, and D, H, and CL fitted well for the subtotals and total biomass for Mongolian oak. Further, AMS0 performed slightly better than AMS1 and AMS3 for most diameter classes, except for the largest diameter classes (D > 30 cm) in which AMS3 appears to give better predictions.  the other two aggregative biomass modeling approaches. Hence, the MPE, MAE, and MAPE across diameter classes were visualized in order to enhance the quality of validating tree biomass models ( Figure 4). Results show that the three aggregative systems based on D, D and H, and D, H, and CL fitted well for the subtotals and total biomass for Mongolian oak. Further, AMS0 performed slightly better than AMS1 and AMS3 for most diameter classes, except for the largest diameter classes (D > 30 cm) in which AMS3 appears to give better predictions.

Comparison of Biomass Models
In the present study, we compared our locally developed equations both based on (1) D, and (2) D and H with the previously developed allometry formulas for Mongolian oak from China State Forestry Administration (SFA [42]) and Wang [36] (Figure 5). The relationship between observed and predicted biomass was analyzed for all categories (i.e., root, stem, branch, foliage, crown, aboveground, and total biomass). Both of the SFA [42] and Wang [36] biomass models were correlated with higher prediction errors compared with the models developed in the present study, verified by quantified measures of relative error (RE) and root mean square error (RMSE). The SFA's equations hugely overestimated branch, foliage, and crown equations by as much as 73.13%, 36.37%, and 58.91% for the D models, and 65.73%, 28.35%, and 51.97% for the D and H models, respectively, while Wang's [36] equations gave a notable underestimation on the root, branch, and crown biomass-as much as 7.14%, 12.50%, and 12.84% (Figure 5), respectively. More importantly, the present study yielded smaller RMSE values in most categories for both D and D and H biomass equations than did the other two previous studies. Finally, some notably large RMSE values occurred in four categories for the D models of the earlier studies: namely RMSE for SFA [42]

Comparison of Tree Carbon Quantification Methods
The RE was computed for comparing the quantification of tree components, subtotals, and total carbon stocks using the six different procedures ( Figure 6). Apart from foliage, the result indicates that there were three distinct strata of RE formed in each category (component, subtotal and total). RRCCc and RWMCc appeared to be the best methods with the smallest error and deserved to be called as the first stratum. Subsequently, three methods followed as the second stratum, namely, RCMCc, TWMCc, and GCCCf-1. GCCCf-2 was found as the third stratum for having the biggest error and undeniably appearing as the worst method. The average RE of all categories combined for the first, second, and third stratum was 3.01 ± 0.18%, 4.08 ± 0.09%, and 10.57%, respectively. Overall, the six approaches for quantifying carbon stock can be ranked based on the relative error value as follows: RRCCc < RWMCc < RCMCc < TWMCc < GCCCf-1 < GCCCf-2.

Comparison of Tree Carbon Quantification Methods
The RE was computed for comparing the quantification of tree components, subtotals, and total carbon stocks using the six different procedures ( Figure 6). Apart from foliage, the result indicates that there were three distinct strata of RE formed in each category (component, subtotal and total). RRCCc and RWMCc appeared to be the best methods with the smallest error and deserved to be called as the first stratum. Subsequently, three methods followed as the second stratum, namely, RCMCc, TWMCc, and GCCCf-1. GCCCf-2 was found as the third stratum for having the biggest error and undeniably appearing as the worst method. The average RE of all categories combined for the first, second, and third stratum was 3.01 ± 0.18%, 4.08 ± 0.09%, and 10.57%, respectively. Overall, the six approaches for quantifying carbon stock can be ranked based on the relative error value as follows: RRCCc < RWMCc < RCMCc < TWMCc < GCCCf-1 < GCCCf-2.

Discussion
In the present study, the development of biomass allometric models and the analysis of carbon concentrations were based on a set of data having a broad range of tree diameters and a wide geographic coverage, allowing estimations to be more precise. Specifically, this comprehensive sampling strategy enabled the effects of various tree components on carbon concentrations of Q. mongolica in five different growing regions, along with their variation within tree sizes, to be thoroughly and duly scrutinized. Carbon concentrations were found to be significantly different within tree components and tree-growing regions, ranging from 0.4194 ± 0.0171 in the root of ZGCME region to 0.4858 ± 0.0177 in the foliage of XXMS region. This finding is in accordance with those that have been reported elsewhere in the literature, which show significant differences in carbon concentration in components of the balsam fir in Canada [43], 10 Chinese temperate tree species in northeast China [37], and eight subtropical tree species in southern China [44]. Moreover, Rodríguez-Soalleiro et al. [34] confirmed the importance of biomass components and tree species as the source of variation in carbon concentration of several broadleaved species, which were found to explain about 52% and 17% of the variance, respectively. Overall, carbon concentration was found highest in foliage (0.4670), lowest in root (0.4406), and slightly lower in branches (0.4491) than stem (0.4568). Thomas and Malczewski [26] and Zhang et al. [37] studied the carbon content of the same specific species in northeast China. Thomas and Malczewski [26] focused on analyzing the stem wood carbon content, while Zhang et al. [37] reported the carbon content in all tree components. By comparison, the total range of carbon concentration across whole components in this study was figured to be more stable and narrower than reported by Zhang et al. [37] (0.4620-0.5140), and the stem carbon concentration in this study was about 4% smaller than Thomas and Malczewski [26].
To date, there are two completely opposing ways to calculate carbon stocks. Some researchers have used a constant generic 0.50 carbon concentration value to be multiplied by a biomass value [45][46][47], in line with the recommendation of Roy et al. [48] and Wenzl [49]. On the contrary, some

Discussion
In the present study, the development of biomass allometric models and the analysis of carbon concentrations were based on a set of data having a broad range of tree diameters and a wide geographic coverage, allowing estimations to be more precise. Specifically, this comprehensive sampling strategy enabled the effects of various tree components on carbon concentrations of Q. mongolica in five different growing regions, along with their variation within tree sizes, to be thoroughly and duly scrutinized. Carbon concentrations were found to be significantly different within tree components and tree-growing regions, ranging from 0.4194 ± 0.0171 in the root of ZGCME region to 0.4858 ± 0.0177 in the foliage of XXMS region. This finding is in accordance with those that have been reported elsewhere in the literature, which show significant differences in carbon concentration in components of the balsam fir in Canada [43], 10 Chinese temperate tree species in northeast China [37], and eight subtropical tree species in southern China [44]. Moreover, Rodríguez-Soalleiro et al. [34] confirmed the importance of biomass components and tree species as the source of variation in carbon concentration of several broadleaved species, which were found to explain about 52% and 17% of the variance, respectively. Overall, carbon concentration was found highest in foliage (0.4670), lowest in root (0.4406), and slightly lower in branches (0.4491) than stem (0.4568). Thomas and Malczewski [26] and Zhang et al. [37] studied the carbon content of the same specific species in northeast China. Thomas and Malczewski [26] focused on analyzing the stem wood carbon content, while Zhang et al. [37] reported the carbon content in all tree components. By comparison, the total range of carbon concentration across whole components in this study was figured to be more stable and narrower than reported by Zhang et al. [37] (0.4620-0.5140), and the stem carbon concentration in this study was about 4% smaller than Thomas and Malczewski [26].
To date, there are two completely opposing ways to calculate carbon stocks. Some researchers have used a constant generic 0.50 carbon concentration value to be multiplied by a biomass value [45][46][47], in line with the recommendation of Roy et al. [48] and Wenzl [49]. On the contrary, some authors prefer to multiply the biomass value with the total weighted mean carbon concentration (TWMCc) of an individual tree [50,51], while others prefer to used respective mean component carbon concentration (RMCCc) as a multiplying factor to obtain carbon stock estimates [52]. These second forms of approaches are considered a better approach since many studies have reported that carbon concentration varies greatly depending upon the type of tree species and tissues [37,53]. Accordingly, using a generic carbon concentration conversion factor (0.50) may cause serious errors in estimating carbon stock. In the current study, six different alternative procedures of quantifying carbon stock were tested for Mongolian oak. The regional respective component carbon concentration (RRCCc) was the best procedure to estimate the carbon stock of an individual tree though implementing the RRCCc requires the comprehensive data of carbon concentration according to the tissue types and specific growing region. Hence, we propose the usage of a regional weighted mean carbon concentration (RWMCc) procedure as an appropriate alternative approach to quantify carbon stock, which performed equally as well as RRCCc. Our result evidence that employing RWMCc will only give 0.25% higher relative error compared with RRCCc. Our data and analyses firmly recommend against utilizing the generic carbon concentration proportion of 0.50 for Mongolian oak as it will likely lead to significant biases in yielding carbon stock estimation, just as it has been reported in other studies for other tree species [28,32,39].
Tree diameter at breast height (D) is an indispensable predictor in biomass modeling and indubitably appears as among the most fundamental elements in forest inventory practice [54,55]. As a sole predictor in estimating tree biomass, D has always been providing an accurate model prediction for many species across numerous regions and delivering a simple model structure that is acquainted with being easy-to-use for forest practitioners [56,57]. The insertion of tree height (H) and crown characteristics (i.e., crown length [CL] and crown width [CW]) into biomass allometric models has been shown to generate a significant improvement in the model fitting and performance. In the current study, D, D and H, and D, H, and CL were used as a predictor(s) to develop each of the three systems of additive biomass model (Equations (1)-(4)). One of our available crown attribute predictors, CW, was intentionally eliminated in the model development since it was not significantly improving the model fitting performance for all tree component biomass equations. Our results noted that the inclusion of H and H and CL in the stem and foliage biomass equations, respectively, could improve the overall model fitting performance. This finding is in accordance with previous literature [6,15,58,59].
In this study, the aboveground (stem, branch, and foliage) and belowground (root) biomass data were obtained across a range of diameter (1.70-31.90 cm) and height (2.40-27.00 m). Previous research revealed that Mongolian oak might grow to exceed the diameter and height range investigated in this study [60]. However, the existences of those large sizes of Mongolian oak are extremely rare in practice, especially in natural secondary forest. Moreover, the Chinese Government has strictly forbidden to cut the large size trees (D > 50 cm) for any kind of purposes. Thus, our newly developed aggregated biomass model systems are inappropriate for estimating the tree components, subtotals, and total biomass exceeding the diameter and height range provided in the present study.
The three aggregative modeling approaches were constructed to develop biomass models. The AMS0 contains no parameter restriction, while the AMS1 and AMS3 contain one and three parameter restrictions, respectively, which guarantees that the sum of the component estimations will be exactly the same as the direct estimation of subtotals and total biomass. In this study, the focus was on demonstrating the differences between AMS0, AMS1, and AMS3 due to the parameter restrictions imposed on the model systems. Our results indicate that all the three aggregative modeling approaches performed well in the context of average prediction errors for components, subtotals, and total biomass estimations. Compared with AMS1 and AMS3, AMS0 was actually better to predict the component, crown, aboveground, and total biomass, although AMS1 only utilizes the total biomass model as a dependent equation, while AMS3 specifies all of the crown, aboveground, and total biomass. These results conform with the previous findings reported by Zhao et al. [21] and Dong et al. [24], which specifically stated that AMS0 is the most suitable approach to construct aggregated systems of biomass models.
It was realized that a limited number of individual tree biomass equations exist for the particular tree species, especially in Northeast China. Wang [36] developed biomass models for this species from the Maoershan Ecosystem Research Station in Heilongjiang Province, Northeastern China. The biomass models were also fitted by using a range of tree sizes, but unfortunately, the equations were lacking an additivity property for component, subtotal, and total biomass. Compared with the present research, Wang's [36] study produced similar predictions in estimating stem, aboveground, and total biomass for Mongolian oak, while there were some notable differences in estimating root, branch, and crown biomass. Moreover, the present study was compared with the previous equations developed by China State Forestry Administration [42]. The State's equations tend to give a distinguished underestimate of stem biomass and produce a remarkable overestimate prediction on branch and crown equations, chiefly for large-sized trees, whereas delivered similar estimations in estimating root, aboveground, and total biomass. The potential explanations of these discrepancies may be due to the differences in sampling sizes, sites, and the characteristic of stands where the trees were sampled. These issues may inevitably relate to the disparities in terms of climate variations, altitude, slope, and soil conditions as external factors; along with stand density, competition from the neighboring trees, and the morphological traits of tree roots as internal factors [6,61,62].

Conclusions
We developed three aggregative systems of biomass equation for Mongolian oak occupying five mountainous regions in Northeast China: one based on D only, one on the combination of D and H, and one on the combination of D, H, and CL. We found that the inclusion of H into the aggregative systems improved the performances and model fit significantly for the stem, aboveground, and total biomass models, while adding CL into the aggregative model system only gave a slight improvement for the model fitting of the foliage and crown biomass models. All of the new systems developed here guarantee the additivity property between the components, subtotals, and total biomass, and provide more accurate estimations than previous research works. There was not available any single function or structural system that was perfect for estimating biomass for all of the tree components, crown, aboveground, or total biomass. Based on numerous fit statistics, the overall performances of the three aggregative model systems can be ranked in the following order: AMS0 > AMS1 > AMS3.
The carbon concentrations were found to be significantly different among tree components and growing regions in our natural secondary forest sites. Different procedures for quantifying the carbon stocks of Mongolian oak individual trees were evaluated. The results indicate that considering carbon concentration differences among tree tissues and growing regions appeared to be the best approach. Contrarily, the usage of a 0.50 generic carbon concentration value will give a relatively high bias in quantifying carbon stocks and evidently appeared as the worst method. The tree species used in this research is widespread in northeastern China. The aggregative biomass models and the carbon concentration analyses presented in this study might become an essential reference for researchers, forest practitioners, and governments, particularly for the Chinese National Forest Inventory. However, the newly developed model systems are not recommended to be used in estimating the biomass of trees outside the scope of the data, both geographically and in terms of tree size.