Stratification-Based Forest Aboveground Biomass Estimation in a Subtropical Region Using Airborne Lidar Data

Species-rich subtropical forests have high carbon sequestration capacity and play important roles in regional and global carbon regulation and climate changes. A timely investigation of the spatial distribution characteristics of subtropical forest aboveground biomass (AGB) is essential to assess forest carbon stocks. Lidar (light detection and ranging) is regarded as the most reliable data source for accurate estimation of forest AGB. However, previous studies that have used lidar data have often beenbased on a single model developed from the relationships between lidar-derived variables and AGB, ignoring the variability of this relationship in different forest types. Although stratification of forest types has been proven to be effective for improving AGB estimation, how to stratify forest types and how many strata to use are still unclear. This research aims to improve forest AGB estimation through exploring suitable stratification approaches based on lidar and field survey data. Different stratification schemes including non-stratification and stratifications based on forest types and forest stand structures were examined. The AGB estimation models were developed using linear regression (LR) and random forest (RF) approaches. The results indicate the following: (1) Proper stratifications improved AGB estimation and reduced the effect of underand overestimation problems; (2) the finer forest type strata generated higher accuracy of AGB estimation but required many more sample plots, which were often unavailable; (3) AGB estimation based on stratification of forest stand structures was similar to that based on five forest types, implying that proper stratification reduces the number of sample plots needed; (4) the optimal AGB estimation model and stratification scheme varied, depending on forest types; and (5) the RF algorithm provided better AGB estimation for non-stratification than the LR algorithm, but the LR approach provided better estimation with stratification. Results from this research provide new insights on how to properly conduct forest stratification for AGB estimation modeling, which is especially valuable in tropical and subtropical regions with complex forest types.


Introduction
The Chinese subtropical regions with rich tree species have a higher carbon sequestration capacity than tropical and temporal forests in the rest of Asia and regions at the same latitude in Europe, Africa, and North America [1] and play important roles in regional and global carbon regulations and climate changes [2,3]. These subtropical regions have complex topographies with mountains, hills, and plains, as well as high forest coverage but fragmental and diverse forest patches due to intense disturbances from humans and nature [4,5]. In particular, a large proportion of young plantations with relatively low biomass density but high growth rates have been extensively distributed, becoming an important research hotspot for carbon cycling [5][6][7][8][9]. Therefore, there is a need to understand the spatial patterns of forest biomass distribution in a timely manner.
The advantages of remote sensing technologies in data collection, coverage, and representation make it an important tool for forest aboveground biomass (AGB) estimation [10]. A large number of studies have been conducted using different sensor data such as optical sensors, radar, and lidar (light detection and ranging) (see review papers [10][11][12][13][14]). Optical sensor data, such as Landsat, have been extensively applied to different applications such as forest mapping and AGB estimation [15][16][17], but the data saturation problem, such as greater than 150 Mg/ha, makes AGB estimation difficult for forest sites with high AGB, [18][19][20]. Although a radar data system penetrates a forest canopy to a certain degree depending on the wavelengths, its serious effects from topography, its poor ability in forest classification, as well as the data saturation problem make it unsuitable for forest AGB estimation, especially in mountainous regions [12,21,22]. The Lidar data method with its ability to capture forest vertical features, that is, tree height, effectively solves the data saturation problem in optical and radar data and provides better AGB estimation than other sensor data systems [23][24][25][26].
The studies on lidar-based AGB estimation are mainly based on the following two broad aspects: selection of proper variables from lidar data and establishment of AGB estimation models [10,14]. The common variables from lidar point clouds or canopy height model (CHM) data include mean height, standard deviation, and height percentile (10th to 100th) [27][28][29]. However, not all variables are needed because of their high correlations to each other. It is necessary to identify the key variables for developing AGB estimation models. In general, stepwise regression is often used to identify key variables when a linear regression (LR) approach is used [30]. Sometimes, the relationship between AGB and variables are nonlinear. In this case, random forest (RF) approach provides the ranking of variable importance [31,32]. In reality, there is no universal conclusion with respect to which kinds of variables perform better; it depends on the characteristics of the study areas and specific tree species or forest types [29,[33][34][35]. In addition to height-related variables, canopy cover, density, and volume, which enhance the three-dimensional (3D) features, are also used for AGB estimation modeling [27,36].
Forest AGB can be estimated using linear or nonlinear models, or machine learning approaches [10,25]. The LR model is the most frequently used approach, although some previous research also did nonlinear transform using logarithms and square roots for response or explanatory variables [24,37,38]. In reality, the linear-based models are not optimal for AGB estimation because AGB is not linearly related to variables, considering the complexity of tree species composition and the impact of environmental conditions (e.g., topography, soil, and moisture) on tree growth [10,39]. AGB is also related to tree density, canopy structure, and tree growth rate, in addition to tree height [25,32]. In this case, nonparametric algorithms such as support vector machine and neural network provide better estimations [9,10]. In recent years, deep learning, due to its powerful data mining ability, has been employed in AGB estimation [40]. To date, it is still unclear which algorithm provides better AGB estimation. The machine learning algorithms cannot guarantee better prediction results than the traditional LR approach [9].
In order to improve AGB estimation, many efforts have been conducted such as a combination of lidar and optical sensor data, based on the assumption that lidar data provide vertical features and optical sensor data provide horizontal features [24,41]. However, controversial conclusions have been obtained, depending on the characteristics of the study areas and datasets used [24,[42][43][44].

of 22
A similar situation is the combination of lidar and radar data [45,46], which may or may not improve estimation [29,47].
The AGB estimation models based on individual forest types have been proven to be effective for improving estimation accuracy [9,20]. Previous research indicated that the relationships between forest parameters and lidar-derived variables varied, depending on the study areas under investigation, tree species composition, and physical settings [48,49]. The AGB of a single tree is largely determined by the diameter at breast height (DBH), the height, and the wood density, which are closely related to unique properties of the tree species [25]. This provides rationale for AGB estimation by stratifying forest types into different stratums and developing models separately [24,50]. In addition to the stratification of forest types, other stratification approaches have been based on topography and the cosine of sun incident [20,51].
Although the stratification-based AGB estimation approach has been proven to be effective for improving AGB modeling results, its real applications have some limitations because using stratification implies using a lower number of sample plots, and as a result it is difficult to develop AGB models for different strata. In previous research, stratification of forest types was mainly based on the available sample plots and the analyst's experience. It is unclear how many strata can be used and how to separate the forest groups for AGB modeling. Therefore, this research attempts to identify a suitable approach to conduct the stratification of forest types and to identify a proper modeling approach for the selected stratum through comparative analysis of modeling results using LR and RF approaches. It is expected to provide new insights on how to implement forest stratification and how to select suitable modeling approaches corresponding to specific forest types in a subtropical forest ecosystem, and therefore provide the best AGB modeling results.

Materials and Methods
We selected the Gaofeng Forest Farm as a case study to explore whether different stratification scenarios of forest types improve AGB estimation and which algorithm, i.e., LR or RF, performs better corresponding to a specific forest type. Figure 1 illustrates the framework of this research. It includes the following major steps: (1) collection and organization of AGB samples, (2) calculation of CHM from airborne lidar data, (3) mapping of forest distribution and design of stratification scenarios, (4) identification of key variables and establishment of AGB estimation models, and (5) evaluation and comparison of the modeling results and application of the developed models to the entire study area.

Study Area
Established in 1953, the Gaofeng Forest Farm is located in Nanning, Guangxi Zhuang Autonomous Region and is the largest state-owned farm in Guangxi. The farm sits on the Daming Mountain, which is characterized by hilly terrain that varies from high elevation in the northeast to low elevation in the southwest. This region has a humid subtropical monsoon climate with adequate daylight and abundant year-round precipitation. The average annual temperature is about 21 • C, and the annual average precipitation is 1200 to 1500 mm. The dominant soil in this area is lateritic red soil with deep soil layers. The excellent natural conditions facilitate the growth of various tropical and subtropical species.
The study area is located in Jiepai and Dongsheng within the Gaofeng Forest Farm ( Figure 2) with an area of 56.3 km 2 . According to the lidar-derived digital elevation model (DEM) data, the elevation ranges from 73.4 to 463.6 m and the slope ranges from 0 • to 69.7 • with the majority of slopes between 20 • and 35 • . According to the forest inventory data of 2016, the main forest types in the area are eucalyptus, Chinese fir, Masson pine, and star anise. In particular, eucalyptus plantations have expanded rapidly since 2002 because of their important economic value for local government and farmers. We selected the Gaofeng Forest Farm as a case study to explore whether different stratification scenarios of forest types improve AGB estimation and which algorithm, i.e., LR or RF, performs better corresponding to a specific forest type. Figure 1 illustrates the framework of this research. It includes the following major steps: (1) collection and organization of AGB samples, (2) calculation of CHM from airborne lidar data, (3) mapping of forest distribution and design of stratification scenarios, (4) identification of key variables and establishment of AGB estimation models, and (5) evaluation and comparison of the modeling results and application of the developed models to the entire study area.

Field Survey and Biomass Calculation at the Plot Level
A field survey of the Gaofeng Forest Farm was conducted in January and February 2018, the same period as the airborne lidar data collection. Using the forest distribution mapped in 2016, sample plots were allocated using a stratification sampling approach according to areas of specific forest types and their age groups. From our study area which consisted of two subfarms (Jiepai and Dongsheng) within this forest farm (see Figure 2), a total of 71 sample plots were allocated. After carefully examining the location and measurement quality, 11 samples were eliminated because the forests, mainly eucalyptus due to its short harvest rotation, had been clear cut. Thus, 60 sample plots, including 37 eucalyptus, 17 coniferous (mainly Chinese fir), and 6 other broadleaf forests, were used. Each sample plot contained the same forest type without a mix of different land covers, and therefore each sample represented a specific forest type. Different plot sizes of 20 x 20 m and 25 x 25 m were used during the field survey based on forest ages and working loads. The coordinates of each sample plot were recorded. Tree species, height, crown diameter and height, and DBH for each individual tree within a sample plot were recorded and measured. The allometric models for major tree species in Guangxi (eucalyptus, Masson pine, Chinese fir, hardwood broadleaf species, and softwood broadleaf species) [52] were used to calculate individual tree AGB. The AGB of a sample plot was the total of individual tree AGBs within the sample plot. Then, the AGB, at the plot level, was converted to AGB at one hectare (Mg/ha). The allometric model can be expressed as Remote Sens. 2020, 12, 1101

of 22
where W is the AGB of an individual tree (kg); D is DBH (cm); H is tree height (m); and c1, c2, and c3 are coefficients for a specific tree species. Table 1 lists the coefficients used in Equation (1) for different tree species [52].

Field Survey and Biomass Calculation at the Plot Level
A field survey of the Gaofeng Forest Farm was conducted in January and February 2018, the same period as the airborne lidar data collection. Using the forest distribution mapped in 2016, sample plots were allocated using a stratification sampling approach according to areas of specific forest types and their age groups. From our study area which consisted of two subfarms (Jiepai and Dongsheng) within this forest farm (see Figure 2), a total of 71 sample plots were allocated. After carefully examining the location and measurement quality, 11 samples were eliminated because the forests, mainly eucalyptus due to its short harvest rotation, had been clear cut. Thus, 60 sample plots, Our collaborator conducted the fieldwork from October to November 2017 and provided the coordinates and forest volume, but not the AGB, for each sample plot. A total of 98 sample plots (plot size of 30 x 30 m) were collected. By examining these sample plots and our lidar data, three plots were excluded because they were outside our study area. Thus, the 95 sample plots included 18 eucalyptus, 20 Chinese fir, 29 Masson pine, and 28 other broadleaf forests.
In order to use the field survey data, it was necessary to convert forest volume to AGB for each sample plot. Thus, based on the sample plots from 2018, we used the same approach to calculate tree volume using the DBH and height for each tree species, and then calculated forest volume for each plot by taking the total volume of individual trees within the same plot. On the basis of the relationship between the AGB and the forest volume at plot level, a LR equation for each forest type was developed, in which AGB is a dependent variable and forest volume is an independent variable. The results Remote Sens. 2020, 12, 1101 6 of 22 showed that the R 2 (coefficient of determination) and relative root mean square error (RMSEr) were, respectively, 0.99 and 3.2% for eucalyptus, 0.83 and 7.9% for coniferous species, and 0.92 and 8.4% for other broadleaf species. Thus, all sample plots collected in 2017 were converted from volume to AGB. Because the gap between the two surveys was only two months, we combined the surveys into one dataset to enlarge our sample size. Table 2 summarizes the statistics of samples, showing that the AGB values for all forest types span from 14.7 to 318.6 t/ha with a coefficient of variation (CV) of 0.51. The detailed AGB distribution of different forest types is illustrated in Figure 3. Eucalyptus and Chinese fir have the most sample plots; however, the AGB values are very concentrated within the range of 50 to 00 t/ha, in particular, Chinese fir has the smallest AGB range and standard deviation. Overall, Chinese fir has the smallest CV value of 0.19 and star anise has the largest CV value of 0.57.    In order to develop AGB estimation models for stratified forest types, an accurate forest 7 distribution map is needed. In this study, the Gaofeng Forest Farm provided a forest distribution 8 map that had been developed in 2016 based on a field survey. According to this map and our research

Update of Forest Distribution Map
In order to develop AGB estimation models for stratified forest types, an accurate forest distribution map is needed. In this study, the Gaofeng Forest Farm provided a forest distribution map that had been developed in 2016 based on a field survey. According to this map and our research objective, a classification system of 10 land covers was determined as follows: Chinese fir, Masson pine, eucalyptus, star anise, Mytilaria laosensis, Castanopsis hystrix, Erythrophleum fordii, other broadleaf tree species, bamboo, and non-forest types (e.g., shrubs, impervious surface, water, and bare soil). The forest distribution map with digital vector format was superimposed on the high-resolution aerial photograph (0.2 × 0.2 m) taken at the same time as the airborne Lidar data were gathered, in January 2018. The analyst carefully checked each forest patch by visual interpretation and comparison between the forest distribution map and aerial photograph. If the boundary of a forest patch had changed, a new boundary was drawn, and the corresponding forest type was updated. This visual interpretation ensured high accuracy in updating the forest distribution map and was effective when there was little change in the land covers. Figure 2c illustrates the spatial distribution of the updated forest types that was used for AGB estimation in this research.

Design of Stratification Scenarios
Most previous studies on forest AGB estimation have not differentiated forest types in developing estimation models, ignoring the variations in the relationships between AGB and canopy height among forest types or tree species, largely due to the limitation in collecting a sufficient number of sample plots [50]. A few studies have proven that separately developing AGB estimation models based on stratified forest types improved predictions [24,50]. However, how to stratify forest types, that is, what criteria can be based on and how many strata are optimal, remains unclear. Thus, we designed different stratification scenarios (Table 3) to examine how they affect the estimation results. The first stratification approach was species based, and forest samples were grouped into two broad groups (coniferous and broadleaf) and five finer groups (Masson pine, Chinese fir, eucalyptus, star anise, and other broadleaf species). The second stratification approach was based on the relationships between canopy height and forest AGB ( Figure 4). The red dashed line in Figure 4 divides all forest samples into two groups. Group 1 mainly consisted of eucalyptus, Chinese fir, Castanopsis hystrix, Erythrophleum fordii, and Mytilaria laosensis, and Group 2 mainly consisted of Masson pine, star anise, and other broadleaf species. Although both Chinese fir and Masson pine are coniferous, the relationships between canopy height and AGB were significantly different; in contrast, Chinese fir and eucalyptus showed similar trends in their relationships between canopy height and AGB. As a comparison, all forest types were used as one population for the AGB modeling. Table 3. Different scenarios for the biomass modeling.

Scenarios Description
Non-stratification (NonS) All forest types were used as one population Stratification based on forest types and tree species (SBFT) (1) Two forest groups (SBFT2), coniferous and broadleaf forests (2) Five forest groups (SBFT5), Masson pine, Chinese fir, eucalyptus, star anise, and other broadleaf tree species Stratification based on forest stand structure (SBFSS) Two groups were selected, based on their relationships between canopy height and aboveground biomass (see Figure 4) Remote Sens. 2020, 12, 1101 8 of 22 Stratification based on forest types and tree species (SBFT) (1) Two forest groups (SBFT2), coniferous and broadleaf forests (2) Five forest groups (SBFT5), Masson pine, Chinese fir, eucalyptus, star anise, and other broadleaf tree species Stratification based on forest stand structure (SBFSS) Two groups were selected, based on their relationships between canopy height and aboveground biomass (see Figure 4)

Collection and Processing of Airborne Lidar Data
The airborne lidar data for the Gaofeng Forest Farm were acquired on 13 and 30 January 2018, using a Riegl LMS-Q680i lidar scanner mounted on the LiCHyAOS system. The specifications of the lidar sensor are shown in Table 4. The lidar data provider preprocessed the raw data and labeled each echo signal as ground return or nonground return before delivering to the end users. The nonground return points were used to generate the digital surface model (DSM) and the ground points were used for the digital elevation model (DEM) at a pixel size of 1 x 1 m using a triangulated irregular network (TIN) interpolation algorithm. The canopy height model (CHM) data were produced from the differences between the DSM and DEM. This initial CHM recorded the maximum object height within a pixel. It could have been contaminated by artificial objects such as power poles, power lines, or towers. Thus, pixels with values greater than 50 m were checked on the aerial photograph with spatial resolution of 0.2 m. If these pixels were found to have abnormal values, they were replaced with the mean value using a mean filter operation.
A number of predictive variables for the AGB estimation were extracted from the CHM data with a uniform window size of 20 × 20 m. The central point of each plot was used to link the sample plot to the CHM image to extract variables according to the window size. The commonly used statistical metrics such as mean, quadratic mean, and height percentiles (10 to 100) [28,33,36] were extracted. Meanwhile, we also explored texture measures using the gray level co-occurrence matrix (GLCM) within the size of the sample plot (20 × 20 m) to extract textural images from the CHM data. The texture measures included correlation, contrast, dissimilarity, entropy, homogeneity, second-order moment, and variance [53,54]. These metrics and the corresponding forest AGB of sample plots were exported to SPSS for further analysis.

Development of Biomass Estimation Models
Many previous studies have indicated that when lidar data are used for forest biomass estimation, linear regression (LR) is one of the most widely used parametric modeling approaches [34,36,55]. The LR approach assumes a linear relationship between the dependent variable (forest AGB here) and predictive variables (lidar-derived variables here). Because a large number of potential variables can be used for AGB modeling, it is necessary to identify key variables that largely explain the variance of the dependent variable [1]. One common approach to identify key variables is stepwise regression analysis.
Stepwise regression analysis is able to determine inclusion or exclusion of variables based on the test statistics of estimated coefficients via a series of F-tests or T-tests. It builds a model by successively adding or removing one variable at a time [30]. The key variables that were finally selected were used to develop the AGB model using the LR approach. In order to understand the importance of different variables in the developed LR model for AGB estimation, standard coefficients (beta) were used to show which variables were more important than others; that is, the higher beta value (absolute value) indicated its higher importance in the LR model.
Considering the potentially nonlinear relationships between AGB and lidar-derived metrics, RF is an alternative. Approach for identifying key variables for AGB modeling. RF is an improved decision tree-based algorithm that uses the bagging and boosting theory. Many trees are produced from the random selection of variables and samples, and thus the final result is based on the ensemble of the individual trees [31]. Compared with the traditional decision tree algorithm, the RF algorithm has some advantages such as insensitivity to noisy data in training datasets, using discrete or continuous datasets, and dealing with large datasets [56]. Many previous studies have indicated that the RF algorithm provides better classification results or estimations than the decision tree algorithm; therefore, the RF algorithm has been widely applied to quantitative analysis such as land-cover classification and AGB estimation [9,15,57,58]. We used the RF algorithm to identify potential key variables, then correlation analysis was used to examine the relationships between these potential variables. For a variable having high correlation with other variables while having less importance in the ranking, this variable was removed, and the RF algorithm was re-run. This procedure was repeated until the minimum number of variables was obtained without reducing the AGB modeling performance. The selected key variables were finally used to develop the AGB estimation model using a RF approach.

Evaluation of Biomass Modeling Results and Application of the Developed Models to Entire Study Area
An evaluation of the modeling performances is often required to understand whether the developed model is suitable for the AGB estimation for the entire study area [10]. Different quantitative measures such as R 2 , RMSE, mean absolute error, system error, Akaike's information criterion, and Bayesian information criterion can be used [59,60]. In addition to the traditional accuracy assessment measures, Valbuena et al. [61] indicated the necessity to also conduct an accuracy assessment of the AGB estimates through hypothesis testing and overfitting evaluation. Considering the validation data used, two methods are often used for the evaluation of AGB estimates [18]. The first approach uses independent testing samples that are different from training samples. The second approach is cross-validation, which is one of the most widely used tools for accuracy assessment, especially with a small number of samples [62]. Cross-validation randomly partitions the samples into several disjointed subsets of approximately equal size (m-fold); each subset is used as testing data, and the remaining subsets are used as training data. The model generated from the training subsets is applied to the testing data to estimate the expected discrepancy. This process is repeated until all possible subsets of samples have been used once; then, the estimation accuracies across all blinded tests are combined to give an overall performance estimate [24,63]. In this study, we used leave-one-out cross-validation for model development and accuracy assessment. One sample was left out and the remaining samples were used to construct the LR model or RF model. This process was repeated, removing each sample in the sample data one by one. Therefore, for each validation, almost all samples were used as training data, which produced reliable prediction, avoiding the impact from random factors. The accuracy of models was assessed by R 2 , RMSE, and RMSEr [64]. By comparing the performances of those models developed using LR and RF algorithms under four stratification scenarios, the best model for the corresponding scenario was selected and applied to predict the forest AGB of the entire study area. In addition to the overall accuracy assessment, R 2 , RMSE, and RMSEr were calculated for evaluation of the modeling results according to individual forest types (i.e., Masson pine, Chinese fir, eucalyptus, star anise, and others) and AGB ranges (i.e., <50, 50-100, 100-150, 150-200, 200-250, and >250 Mg/ha). On the basis of these evaluation results, the best modeling results for the forest types were combined to generate a new AGB map with the highest estimation accuracy.

Comparative Analysis of Biomass Estimation Models
A comparative analysis of the AGB models (see Table 5) indicates that in the LR-based AGB models, (1) H ME is often selected for the AGB estimation models and has a higher beta value than other selected variables, implying that it is a more important variable and (2) for individual tree species such as Chinese fir, eucalyptus, and star anise with relatively simple stand structures, only one variable is selected, but for the complex forest types such as coniferous or broadleaf forests, more variables related to forest stand structure, in addition to H ME , are needed. Table 5 also shows that overall the RF-based models have higher R 2 values and require more variables than LR-based models. The selected variables vary in different stratification scenarios. Even in the same scenario of SBFT5, the selected variables for each forest type vary, implying the importance of identifying suitable variables for specific tree species.

Overall Evaluation of Biomass Modeling Results
The AGB estimation results based on different scenarios (Table 6) indicate that the LR model has slightly poorer performance than the RF model based on NonS and SBFT2, but the inverse based on the SBFT5 and SBFSS scenarios. Overall, the LR model based on the SBFT5 scenario performs the best, with R 2 of 0.8 and RMSE of 25.36 Mg/ha; in contrast, the LR model based on the SBFT2 scenario has the poorest performance, with R 2 of 0.43 and RMSE of 44.07 Mg/ha, implying the importance of detailed stratification of forest types. The results in Table 6 also indicate that improper stratification, i.e., coniferous and broadleaf forest (SBFT2), cannot improve the AGB estimation, but a good stratification, such as the SBFSS scenario here, can significantly improve the AGB estimation, that is, R 2 increased from 0.43 (based on SBFT2) to 0.78 (based on SBFSS) and RMSE decreased from 44.07 Mg/ha to 26.0 Mg/ha. On the one hand, more detailed stratification (SBFT5) performs similarly to the SBFSS scenario when the LR model is used. On the other hand, when the RF model was used, the SBFSS scenario provided the best AGB modeling results. This situation implies that very detailed stratification of forest types is not needed. The finding is valuable for guiding the collection of sample plots and the possibility of using a small number of them for the AGB estimation in some forest types. Note: R 2 , coefficient of determination; RMSE, root mean squared error; RMSEr, relative root mean squared error; NonS, SBFT2, SBFT5, and SBFSS are defined in Table 3.
A comparative analysis of the relationships between AGB estimates and reference data ( Figure 5) indicates that (1) both the LR and RF models based on the SBFT5 or SBFSS scenarios (Figure 5a3,a4,b3,b4) perform better than those based on NonS or SBFT2 (Figure 5a1,a2,b1,b2); (2) compared to the the LR modeling results, the RF modeling results have a common problem, that is, when AGB is less thañ 50 Mg/ha, AGB was overestimated; (3) except for the LR modeling result based on the SBFT5 scenario, underestimation is obvious, especially when AGB is greater than 220 Mg/ha. Figure 5 also shows that the LR model based on the SBFSS scenario (Figure 5a4) is the best when AGB is less than~120 Mg/ha, and the LR model based on the SBFT5 scenario (Figure 5a3) is the best when AGB is greater thañ 270 Mg/ha.

Evaluation of Biomass Estimation Results According to Individual Forest Type
The evaluation results according to different forest types (Table 7) indicate that the modeling based on the SBFT5 or SBFSS scenarios using the LR and RF approaches performed much better than modeling based on the NonS and SBFT2 scenarios. The RMSEr values are 19.7% to 21.7% for Masson pine and 15.9% to 18.9% for Chinese fir based on the SBFT5 and SBFSS scenarios, respectively, as compared with 32.6% to 39.5% and 39.2% to 52.6% based on the NonS and SBFT2 scenarios, respectively. Overall, the LR model based on the SBFT5 scenario provided the best results for eucalyptus and star anise, while the LR model based on the SBFSS scenario provided the best results for Chinese fir (RMSEr 16%), and the RF model based on the SBFSS scenario provided the best result for "others" class (other broadleaf species). This evaluation result implies the importance of proper stratification for improving AGB estimation. For example, if AGB models were separately developed for Masson pine and Chinese fir, the RMSE for Masson pine could be reduced from 55.7 Mg/ha (NonS) and 61.4 Mg/ha (SBFT2) to 31.2 Mg/ha (SBFT5), and for Chinese fir from 32.7 and 43.8 Mg/ha to 15.6 Mg/ha, respectively. Table 7 also indicates that for different forest types, the best modeling approach varies, for example, for the SBFT5 scenario, the RF model provided slightly better modeling results for Masson pine and Chinese fir than for LR, whereas the LR model provided better results for eucalyptus, star anise, and "others." Another finding in Table 7 is that star anise provides much better estimation using the LR model based on the SBFT5 scenario, implying that when the number of sample plots is not large, the LR model is more reliable than the RF model. However, for a complex forest type, the "others" class, in this study, with a limited number of sample plots, the SBFSS scenario with an RF algorithm provides the best estimation, implying the advantage of using the SBFSS scenario.  ; (1, 2, 3, and 4) respectively, represent four stratification scenarios NonS, SBFT2, SBFT5, and SBFSS, as defined in Table 3.  ; (1, 2, 3, and 4) respectively, represent four stratification scenarios NonS, SBFT2, SBFT5, and SBFSS, as defined in Table 3.

Evaluation of Biomass Estimation Results According to Biomass Ranges
The evaluation results according to the AGB ranges (Table 8) indicates that proper stratification of forest types (SBFT5 and SBFSS scenarios) considerably improves the AGB estimation for each AGB range, especially when the AGB is low (e.g., less than 50 Mg/ha) or high (e.g., greater than 200 Mg/ha). When the AGB is less than 50 Mg/ha, the LR model based on the SBFT5 scenario provides the best estimation, followed by the SBFSS scenario with RMSE less than 13.5 Mg/ha. When the AGB is 50 to 150 Mg/ha, the LR model based on the SBFSS scenario has the best estimation with RMSE 17.1 to 24.2 Mg/ha, and when the AGB is greater than 150 Mg/ha, the SBFT5 scenario has the best estimation with RMSE 26.7 to 40.5 Mg/ha. Comparing Table 8 with Figure 5, we find that stratification of forest types using the SBFT5 and SBFSS scenarios considerably reduced the underestimation problem. In particular, when the AGB was greater than 200 Mg/ha, the SBFT5 scenario greatly reduced the underestimation problem.
As shown in Table 6, the SBFT5 and SBFSS scenarios have much better modeling estimates than the NonS and SBFT2 scenarios, no matter which algorithm, LR or RF, is used. However, considering different AGB ranges, the LR and RF algorithms perform differently. For instance, based on the SBFT5 scenario, the LR algorithm performs much better than the RF algorithm when the AGB is less than 50 Mg/ha or more than 200 Mg/ha and slightly better than the RF algorithm when the AGB is 50 to 100 or 150 to 200 Mg/ha, but the LR algorithm performs poorly as compared with the RF algorithm when the AGB is 100 to 150 Mg/ha (RMSEr is 26.1% for the LR algorithm and 23.6% for the RF algorithm). This situation implies that for a finer stratification scenario, the overfitting problem with the RF algorithm is obvious, resulting in an overestimation and underestimation problem, especially when the number of samples is relatively small. For the SBFSS scenario, the LR algorithm performs much better than the RF algorithm when the AGB is less than 50 Mg/ha (RMSEr is 38.1% for the LR algorithm and 80.0% for the RF algorithm) and slightly better than the RF algorithm when the AGB is 50 to 200 Mg/ha, (difference within 0.3% to 2.0% between the LR and RF algorithm results), but the LR algorithm performs poorly as compared with the RF algorithm when the AGB is greater than 200 Mg/ha (difference within 1.8% to 3.4% between the LR and RF results). This situation implies that when the AGB is more than 200 Mg/ha, the AGB and canopy height cannot meet the linear assumption, because the SBFSS scenario is grouped based on the relationship between the canopy height and the AGB (see Figure 4). In this case, the RF algorithm handles the nonlinear relationship better than that of the LR algorithm. The results in Table 8 imply that it is necessary to select different algorithms to develop specific AGB estimation models by taking the stratification and AGB ranges into account.

Comparative Analysis of Biomass Spatial Distributions
The spatial distributions of AGB predictions based on the LR and RF models under different scenarios ( Figure 6) indicate that most of the study area has 50 to 200 Mg/ha of AGB, and a few areas have less than 50 Mg/ha or greater than 200 Mg/ha. The spatial distribution indicates that the LR model prediction results have many more pixels with less than 50 Mg/ha or higher than 200 Mg/ha as compared with the RF model prediction results, implying the inability of the RF model to predict very low or very high AGB amounts. Table 7 indicates that different forest types require specific stratification and modeling approaches; thus, the best modeling result for each forest type provides the best AGB prediction map (Figure 7). A comparison of the spatial patterns in Figures 6 and 7 shows that the number of pixels with very low (less than 50 Mg/ha) or very high (greater than 250 Mg/ha) values is obviously large in Figure 7, implying much improved modeling results. Remote Sens. 2019, 11, x FOR PEER REVIEW 6 of 25 Figure 6. A comparison of biomass prediction results. (a and b) represent linear regression and random forest models, respectively; (1, 2, 3, and 4) respectively, represent the four stratification scenarios NonS, SBFT2, SBFT5, and SBFSS, as defined in Table 3. Table 7 indicates that different forest types require specific stratification and modeling approaches; thus, the best modeling result for each forest type provides the best AGB prediction map (Figure 7). A comparison of the spatial patterns in Figures 6 and 7 shows that the number of pixels  ; (1, 2, 3, and 4) respectively, represent the four stratification scenarios NonS, SBFT2, SBFT5, and SBFSS, as defined in Table 3. with very low (less than 50 Mg/ha) or very high (greater than 250 Mg/ha) values is obviously large in Figure 7, implying much improved modeling results. Figure 7. Spatial distribution of aboveground biomass estimation using the best-performing model for each forest type.

Impact of Stratification on Aboveground Biomass Modeling
Forest AGB is related to tree species composition and tree density in a unit, in addition to wood density, DBH, and tree height of individual tree species. For optical sensor data, data saturation is a big problem that results in poor AGB estimation because optical sensor data mainly capture land surface features without vertical forest stand information [20]. Lidar data measure tree or canopy height, and thus are regarded as the most accurate data source for AGB estimation because of the strong correlation between AGB and canopy or tree height [25]. Some previous studies using Landsat data have proven the effectiveness of stratification of forest types for improving AGB estimation [9,20]. However, airborne lidar data are only available at experimental sites with a limited area; thus, sample plots are often limited [24,65] and non-stratification is common in most previous studies. Cao et al. [27] confirmed that stratification of forest types reduced the overall RMSEr from 22.7% to 15.84% for coniferous forest type, 15.41% for broadleaf forest type, and 18.33% for mixed forest type. Our research indicates that proper stratification of forest types considerably improved the AGB modeling performance, but improper stratification reduced performance. For example, although Masson pine and Chinese fir belong to coniferous forests, their crown size and shape and their relationships between AGB and canopy height are different. This research shows that the AGB modeling based on coniferous forest (e.g., a combination of Masson pine and Chinese fir) reduced modeling performance, while modeling based on individual tree species considerably improved it. As indicated in Table 7, the RMSEr values for Masson pine and Chinese fir using the LR model based on the SBFT2 scenario increased by 3.7% and 13.4%, respectively as compared with those based on the NonS scenario, while they reduced by 24.5% and 17.1% based on the SBFT5 scenario. The RF-based

Impact of Stratification on Aboveground Biomass Modeling
Forest AGB is related to tree species composition and tree density in a unit, in addition to wood density, DBH, and tree height of individual tree species. For optical sensor data, data saturation is a big problem that results in poor AGB estimation because optical sensor data mainly capture land surface features without vertical forest stand information [20]. Lidar data measure tree or canopy height, and thus are regarded as the most accurate data source for AGB estimation because of the strong correlation between AGB and canopy or tree height [25]. Some previous studies using Landsat data have proven the effectiveness of stratification of forest types for improving AGB estimation [9,20]. However, airborne lidar data are only available at experimental sites with a limited area; thus, sample plots are often limited [24,65] and non-stratification is common in most previous studies. Cao et al. [27] confirmed that stratification of forest types reduced the overall RMSEr from 22.7% to 15.84% for coniferous forest type, 15.41% for broadleaf forest type, and 18.33% for mixed forest type. Our research indicates that proper stratification of forest types considerably improved the AGB modeling performance, but improper stratification reduced performance. For example, although Masson pine and Chinese fir belong to coniferous forests, their crown size and shape and their relationships between AGB and canopy height are different. This research shows that the AGB modeling based on coniferous forest (e.g., a combination of Masson pine and Chinese fir) reduced modeling performance, while modeling based on individual tree species considerably improved it. As indicated in Table 7, the RMSEr values for Masson pine and Chinese fir using the LR model based on the SBFT2 scenario increased by 3.7% and 13.4%, respectively as compared with those based on the NonS scenario, while they reduced by 24.5% and 17.1% based on the SBFT5 scenario. The RF-based model based on the NonS and SBFT2 scenarios provided similar modeling performances for Masson pine and Chinese fir, but based on the SBFT5 scenario, the RMSEr values decreased by 13.3% and 22.8%, respectively for those species. This result differs from some previous studies [27,29,39]. The possible reason is that Chinese fir is often in a pure plantation with homogenous stand structure and age, while the stand structure of Masson pine is often complex and mixed with some broadleaf tree species in the canopy and with shrubs and grass under canopy, as shown in Figure 8. Our research indicates that the SBFT5 scenario provides the best modeling performance but requires a large number of sample plots for each tree species. This is often a challenge for most studies, considering the cost and labor of sample data collection. Another challenge in using the SBFT5 scenario is the detailed classification of tree species, which is often difficult because of the similar spectral signatures of different broadleaf tree species. However, many tree species have similar relationships between tree height and AGB, and we can combine them into one group. This implies that we are able to collect fewer sample plots for AGB modeling and still produce similar modeling results, as shown in Table 6, that is, R 2 and RMSE are similar for the SBFSS and SBFT5 scenarios using the LR model, and slightly better when using the RF model. Another advantage is that for some forest types with limited area, the SBFSS scenario can provide modeling results, but the SBFT5 scenario cannot. The challenge of using the SBFSS scenario is to clearly understand the relationships between tree height and AGB for different tree species, and thus expert knowledge is essential for determining the grouping if sufficient sample plots are not available.

Impacts of Sample Sizes on Aboveground Biomass Modeling
The collection of a sufficient number of sample plots is one of the most important steps in AGB estimation [10]. However, in reality, the collection of sample plots is greatly hindered by cost and labor intensity, which seriously affect modeling performance, that is, they affect the determination of stratification of forest types, selection of modeling variables and algorithms, and the evaluation of modeling results. For example, with the numbers of samples for Masson pine, Chinese fir, and eucalyptus at 29, 37, and 55, we developed the LR model for each species. As shown in Table 7, the RMSEr values for Masson pine, Chinese fir, and eucalyptus were, respectively, reduced from 35.8%, 39.2%, and 31.3% based on the NonS, to 20.0%, 18.7%, and 24.8% based on the SBFT5 scenario. However, if the AGB range of a sample plot is not representative, the selection of modeling algorithms is seriously affected. For instance, star anise has only 13 sample plots and only one with the AGB greater than 200 Mg/ha; as a result, the RF model's RMSEr is 56.8% as compared with 21.0% in the LR model. A similar situation occurs for eucalyptus, that is, of 55 samples, only five have AGB greater than 150 Mg/ha. Thus, the RF-based models performed poorly as compared with the LRbased models (Table 7). This implies that collection of a sufficient number of sample plots with good representation of different AGB ranges is critical for AGB estimation.

Impacts of Forest Types on Selection of Modeling Variables and Algorithms
There is no universal conclusion about what variables should be used for AGB modeling, because they depend on forest types, complexity of the landscape under investigation, individual Our research indicates that the SBFT5 scenario provides the best modeling performance but requires a large number of sample plots for each tree species. This is often a challenge for most studies, considering the cost and labor of sample data collection. Another challenge in using the SBFT5 scenario is the detailed classification of tree species, which is often difficult because of the similar spectral signatures of different broadleaf tree species. However, many tree species have similar relationships between tree height and AGB, and we can combine them into one group. This implies that we are able to collect fewer sample plots for AGB modeling and still produce similar modeling results, as shown in Table 6, that is, R 2 and RMSE are similar for the SBFSS and SBFT5 scenarios using the LR model, and slightly better when using the RF model. Another advantage is that for some forest types with limited area, the SBFSS scenario can provide modeling results, but the SBFT5 scenario cannot. The challenge of using the SBFSS scenario is to clearly understand the relationships between tree height and AGB for different tree species, and thus expert knowledge is essential for determining the grouping if sufficient sample plots are not available.

Impacts of Sample Sizes on Aboveground Biomass Modeling
The collection of a sufficient number of sample plots is one of the most important steps in AGB estimation [10]. However, in reality, the collection of sample plots is greatly hindered by cost and labor intensity, which seriously affect modeling performance, that is, they affect the determination of stratification of forest types, selection of modeling variables and algorithms, and the evaluation of modeling results. For example, with the numbers of samples for Masson pine, Chinese fir, and eucalyptus at 29, 37, and 55, we developed the LR model for each species. As shown in Table 7, the RMSEr values for Masson pine, Chinese fir, and eucalyptus were, respectively, reduced from 35.8%, 39.2%, and 31.3% based on the NonS, to 20.0%, 18.7%, and 24.8% based on the SBFT5 scenario. However, if the AGB range of a sample plot is not representative, the selection of modeling algorithms is seriously affected. For instance, star anise has only 13 sample plots and only one with the AGB greater than 200 Mg/ha; as a result, the RF model's RMSEr is 56.8% as compared with 21.0% in the LR model. A similar situation occurs for eucalyptus, that is, of 55 samples, only five have AGB greater than 150 Mg/ha. Thus, the RF-based models performed poorly as compared with the LR-based models ( Table 7). This implies that collection of a sufficient number of sample plots with good representation of different AGB ranges is critical for AGB estimation.

Impacts of Forest Types on Selection of Modeling Variables and Algorithms
There is no universal conclusion about what variables should be used for AGB modeling, because they depend on forest types, complexity of the landscape under investigation, individual tree species, and the algorithms [9,20]. Our research confirms that different variables must be used for specific stratification scenarios and modeling methods. This implies that developing an optimal model should take into account the specific forest types and their complexities of forest stand structure, as well as the modeling algorithm. As shown in Table 5, for a forest composed of a single tree species, such as Chinese fir, eucalyptus, or star anise, where all trees are the same age, only one variable was selected for the LR model, whereas two or more variables were selected for complex forest types such as Masson pine (pine dominant in top canopy and broadleaf tree species dominant under canopy), and "others" (composition of different tree species). For the RF models, more variables were selected than for the LR models. However, using more variables does not guarantee better prediction results because of the overfitting problem in the RF algorithm. As shown in Figure 5, the RF algorithm did not have a good extrapolation ability, especially for the forest types with a limited number of sample plots with low or high AGB ranges. This research shows the importance of selecting proper variables by taking into account the different characteristics of forest types.
The complexity of a forest stand structure also affects the selection of modeling algorithms. As shown in Table 7, the RF algorithm performs better than the LR algorithm in NonS and SBFT2 scenarios, while the inverse occurs in SBFT5 and SBFSS scenarios, implying that for complex forest types, a machine learning algorithm is more suitable than the LR algorithm. This could be due to the fact that the LR algorithm cannot effectively delineate the nonlinear and complex relationships between the forest attributes and the AGB [32,39]. For a forest type with relatively simple stand structure, such as eucalyptus, the LR algorithm performs better, especially for extrapolation when a limited number of samples for high and low AGB samples are used, as shown in Figure 5. This research implies that when the number of sample plots is insufficient, the LR approach is better than the RF approach, and that a sufficient number of samples with relatively low and high AGB must be collected for the RF approach to reduce the over-and underestimation problems [20,24].

Impacts of Amounts of Forest Aboveground Biomass on Modeling Results
How different amounts of forest AGB affect modeling results has not been extensively explored, especially for the Lidar-based AGB estimation. Limited research has examined AGB under-and overestimation using Landsat and ALOS PALSAR data in a subtropical region [9,20,22]. The overestimation in optical or radar data mainly occurs in forest sites with relatively low AGB such as less than 50 Mg/ha, because of the complex composition of land-cover components in one pixel (shrub, grass, or bare soil in addition to tree density), whereas underestimation occurs in dense forests with high AGB, such as more than 150 Mg/ha, due to the data saturation problem [20]. Compared to optical or radar data that mainly capture canopy features without much information on forest height, Lidar data obtain canopy height features that are highly related to AGB, resulting in high estimation accuracy, as shown in this research (Table 8), especially when finer stratification is used. In the SBFSS scenario, when AGB is higher than 200 Mg/ha, the RF approach provides better modeling performance than the LR approach, implying that when AGB is relatively high, AGB may be nonlinearly related to CHM-derived variables. This research also implies that stratification based on forest types and AGB ranges provides better modeling results. This should be explored in the future if a large number of sample plots are available.

Conclusions
This research explored lidar-based AGB modeling in a subtropical forest ecosystem using LR and RF approaches based on different stratifications of forest types. Through comparative analysis of modeling results, this research provided some new insights on how to implement stratification of forest types and how to identify suitable variables for AGB modeling based on different forest types. Specifically, we state the following conclusions: (1) The complex composition and rich tree species in a subtropical ecosystem make AGB estimation a challenge, resulting in high uncertainty in AGB modeling, especially when the AGB is relatively low (e.g., <50 Mg/ha) or high (e.g., >200 Mg/ha). Proper stratification of forest types for AGB modeling is an effective approach to improve AGB estimation. (2) AGB modeling based on specific forest types provides the best AGB modeling results but requires a large number of sample plots for each forest type, resulting in a challenge in sample data collection and in classification of detailed forest types. Stratification based on the relationship between canopy height and AGB (SBFSS in this study area) solves this problem and provides similar modeling performance. (3) Selection of modeling variables depends on the complexity of the forest types or forest stand structures. More variables are often required for the NonS and SBFT2 scenarios or complex forest types (Masson pine) than for detailed stratification with individual forest types. (4) The LR model is suitable for AGB modeling for the forest types with relatively simple forest stand structure such as Chinese fir, eucalyptus, and star anise, while the RF model is suitable for complex forest types such as broad groups, that is, coniferous and broadleaf forest types. (5) Stratification based on the relationship between canopy height and AGB (the SBFSS scenario) is recommended for AGB modeling, considering performance and requirement of sample sizes.