Analysis of Taper Functions for Larix olgensis Using Mixed Models and TLS

Terrestrial laser scanning (TLS) plays a significant role in forest resource investigation, forest parameter inversion and tree 3D model reconstruction. TLS can accurately, quickly and nondestructively obtain 3D structural information of standing trees. TLS data, rather than felled wood data, were used to construct a mixed model of the taper function based on the tree effect, and the TLS data extraction and model prediction effects were evaluated to derive the stem diameter and volume. TLS was applied to a total of 580 trees in the nine larch (Larix olgensis) forest plots, and another 30 were applied to a stem analysis in Mengjiagang. First, the diameter accuracies at different heights of the stem analysis were analyzed from the TLS data. Then, the stem analysis data and TLS data were used to establish the stem taper function and select the optimal basic model to determine a mixed model based on the tree effect. Six basic models were fitted, and the taper equation was comprehensively evaluated by various statistical metrics. Finally, the optimal mixed model of the plot was used to derive stem diameters and trunk volumes. The stem diameter accuracy obtained by TLS was >98%. The taper function fitting results of these data were approximately the same, and the optimal basic model was Kozak (2002)-II. For the tree effect, a6 and a9 were used as the mixed parameters, the mixed model showed the best fit, and the accuracy of the optimal mixed model reached 99.72%.The mixed model accuracy for predicting the tree diameter was between 74.22% and 97.68%, with a volume estimation accuracy of 96.38%. Relative height 70 (RH70) was the optimum height for extraction, and the fitting accuracy of the mixed model was higher than that of the basic model.


Introduction
The taper equation is a mathematical function that studies the diameter of any part of the stem and the corresponding tree height (TH), total TH, and diameter at breast height. The taper function has a wide range of applications in tree stem volume estimation, 3D space model reconstruction, yield estimation, forest planning, and timber production simulation and optimization [1,2]. The taper equation can be used to estimate the diameter at any height on the stem, the height at any diameter, the total stem volume, and the merchantable volume of different specifications [3][4][5]. A simplified Kozak variable-exponential taper function can improve the predictive ability of the model by adding an auxiliary diameter through the mixed effect method [6]. The fixed effect model is simpler than the mixed effect model. However, the mixed model eliminates the error correlation and heterogeneity of hierarchical data, to improve the prediction accuracy of the model and explain the source of random errors [7]. Doyog N D's study of compatible " Max and Burkhart (1976)" showed that the estimated value of the bottom diameter of the trunk is higher than the true value, while the middle and upper diameters are underestimated; the overall height and volume of the stem at a given diameter are overestimated [8]. A compatible system composed of the taper equation, total volume equation and volume sales equation was used for the modeling, and the fitting accuracy of the optimal model reached

Data Source and Processing
In October 2019, nine larch fixed sample plots were selected in the Mengjiagang Forest Farm of Jiamusi city, Heilongjiang Province according to different site qualities, densities, ages and thinning intensities. The DBH, height, CW and relative coordinates of each tree were measured. The DBH of each plot was divided into five grades according to the method of equal cross-sectional area sample trees, and the average DBH of each grade was calculated. Six plots (there were a total of nine plots, but stem analysis was only conducted on six plots due to the weather.) were selected, and five sample trees of different sizes were chosen near the sample plots according to the average DBH of each grade for stem analysis. We collected a total of 30 trees for stem analysis. The DBH, TH, CW, and height of the tree's crown base (HCB) based on the stem analysis were measured. The relative height according to the total tree height was calculated, and the diameter of each relative height was measured, including 0 m, RH2, RH4, RH6, RH8, RH10, RH20, RH25, RH30, RH40, RH50, RH60, RH70, RH75, RH80 and RH90. TLS data were collected by a Trimble TX8 narrow infrared laser beam before stem analysis. Five stations were scanned for each plot, three stations were scanned for each stem analysis, and the scanning time of each station was 3 min. After data preprocessing, a total of 580 larch trees were obtained from 9 plots, and 7887 diameters at different heights were extracted; a total of 451 diameters were extracted from 30 trees for stem analysis. Seventy-five percent of the sample trees were used for model parameterization, and the remainder (25%) were used to evaluate model performance. The parameters of Trimble TX8 are listed in Table 1. A summary of the stem analysis data is provided in Table 2.

Processing Point Cloud Data
The original scanned data were in TZF format. Trimble RealWorks11.1 and LiDAR360 software were used for registration, denoising, 0.01 m minimum point spacing for subsampling, classification of ground points with an improved progressive TIN densification filtering algorithm, establishment of the digital elevation model, normalization by ground class, point cloud segmentation, tree identity matching and other processing. The diameters at different heights were fitted by the least squares method. To reduce the influence of the tree height on the extraction accuracy, the absolute height was replaced by the relative height. The position of stem diameter extraction is shown in Figure 1. When extracting the diameter, the slice thickness was 10 cm. For example, assuming that the diameter at 1.3 m was extracted, point clouds with heights from 1.25-1.35 m were selected for the fitting circle. for the fitting circle.
In this study, the stem analysis data at each relative height of 30 felled trees were selected for verification, which also provided a basis for constructing the taper equation. The analyzed diameter of 30 trees was taken as the reference value, and the corresponding TLS data were taken as the predicted value. The extraction effect was evaluated by the coefficient of determination (R 2 ), root mean square error (RMSE), absolute error (bias) and extraction accuracy (P%). In this study, the stem analysis data at each relative height of 30 felled trees were selected for verification, which also provided a basis for constructing the taper equation. The analyzed diameter of 30 trees was taken as the reference value, and the corresponding TLS data were taken as the predicted value. The extraction effect was evaluated by the coefficient of determination (R 2 ), root mean square error (RMSE), absolute error (bias) and extraction accuracy (P%).

Basic Taper Function
Referring to the latest edition of "Forest Mensuration" (Fourth Edition) [24] and related documents worldwide, six different commonly used types of taper equations were used in this study as candidate models. The relative height values at the lower and upper inflection points of the segmented taper equation were set to RH9 and RH77, respectively [25]. Among the 30 trees in the stem analysis, 75% (23 trees) for model fitting and the remaining 25% (7 trees) for evaluation. TLS data are the same. (1) Simple taper equation: Kozak (1969)-II [26]: Kozak (1969) [26]: Schumacher(1973) [24]: (2) Segmented taper equation: Max and Burkhart (1976) [27]: (3) Variable-exponential taper equation: Weisheng Zeng, Zhiyun Liao (1997) [28]: where d is the stem diameter at a height of h; D is the diameter at breast height; H is the total tree height; h is the measuring point height; a 1 , a 2 , a 3 , · · · , a 9 are the parameters of the model; q = h/H; and t = 1.3/H.

Mixed Effects Models
Because the main purpose of our study was to provide the best taper model to predict the upper stem diameter for point cloud data, the fixed-effects modelling approach guided the model selection procedure. The stem diameter data have a hierarchical structure (several diameter measurements within a single sample tree) and great spatial variability. The taper equation established by the trees in the sample plots can better represent the trends in the stem curve of the forest stand and can also more accurately predict the upper diameter of the trunk and estimate the volume. Once the best model was selected, a nonlinear mixed model with the tree effect was fitted and compared with the fixed-effects model. The mixed model expression is as follows: where Y ij is the dependent variable of the jth observation of the ith tree, which refers to the diameter at the relative height in this paper; m is the number of samples; n i is the number of observations of the ith tree; x ij is the relative height; β is a (p × 1) dimensional vector of the fixed parameters; B ij is a (q × 1) random effects parameter vector; A ij and B ij are the corresponding design matrices; R i is a variance-covariance matrix for the error terms; ε ij is an n-dimensional vector of the residuals; and σ 2 is the variance value.
The following three steps are used to construct the mixed effect model [29]: (1) Determine random parameters. Fit the taper models of different random parameter combinations and compare the fit statistics. Too many parameters may lead to the overparameterization or nonconvergence of the model. Therefore, this study only selects a combination of 1-2 random parameters for fitting and compares the Akaike information criterion (AIC), Bayesian information criterion (BIC) and log-likelihood values of the model fitting. The results show that models with a smaller AIC and BIC generally have a larger logarithmic likelihood value and better fitting quality. (2) Determine the variance-covariance structure within-trees (R). The data of this study has a hierarchical structure among trees, therefore, we needed to solve the problem of the correlation and heteroscedasticity of intra-tree errors. To compensate for the effects of autocorrelation, this study used the first-order autoregressive structure AR (1). This study utilized the work of Davidian (1995), which is common in forestry research, to calculate [30]: where σ 2 refers to the error variance value of the model; Γ i is the time series correlation structure, i.e., the error correlation structure within the tree; and G i is the diagonal matrix of the variance heterogeneity. (3) Determine the structure of the inter-tree variance-covariance (D). The structure of intertree variance-covariance reflects the changes among groups. Taking the Generalized positive-definite matrices commonly used in forestry as an example, the variancecovariance structures of 2 random parameters were selected: where σ 2 u is the variance of the random parameter u, σ 2 v is the variance of the random parameter v, and σ uv is the covariance of random parameters u and v.

Evaluation and Test Models
This study uses independent data that are not involved in modeling for testing. The test of the fixed effects parameter in the mixed effects model is equivalent to the traditional regression analysis test, and the random effect component requires secondary sampling to calculate the random parameter value. Five samples were randomly selected from each tree, and the method of Vonesh was used to calculate the random parameter values [31] whereD is the variance-covariance matrix of the random effect parameters,R i is the variance-covariance structure in the sample tree,Ẑ i is the design matrix, andê i is the actual value minus the predicted value calculated with fixed effect parameters.

Prediction of the Stem Diameter
The taper function can be used to estimate the diameter of any height on the tree. The mixed model was established using 580 trees from 9 plots. It is highly representative and can simulate the stem curve well. It was used in this study to derive the diameter of the upper stem that cannot be obtained by TLS. The prediction of the upper diameter, which was obtained with the mixed model, not only solves the problem of missing point cloud data but also accurately estimates the volume.

Estimated Volume
The ultimate purpose of constructing the taper equation is to calculate the volume. This study adopted the following three types of methods to calculate the volume of 30 trees: Method 1, the optimal taper model with tree effects was integrated numerically to compute the volume; Method 2, the extracted diameters at different relative heights were directly used to calculate the volume with the mid-area sectional measurement method (the diameter that could not be extracted from the upper stem was calculated according to the predicted value by the taper equation); and Method 3, the measured value from the stem analysis at each meter was used to calculate the volume. The measured volume (Method 3) was taken as the reference value and was compared with those of the other two methods to determine their prediction accuracies.

Summary of TLS Extraction Data
The recognition performance of the point cloud data between the ground and the bottom of the tree is limited when establishing the digital elevation model; thus, the difference in the measured diameter of 0 m is also relatively large(to be consistent with the stem analysis, we extracted the point cloud value at 0.1 m as 0 m). Therefore, the extracted diameter at 0 m was compared with the actual measurement. The DBH, TH, HCB, CW and diameters of different positions extracted from 30 trees were analyzed. Table 3 shows a summary of the point cloud-extracted data. When the tree height was between RH75 and RH90, the diameter of some point clouds could not be extracted due to the occlusion caused by branches and leaves. Most of the R 2 values extracted from the diameters at 15 different relative heights were more than 0.92. The extraction results of the DBH and tree height are good, with accuracies of 98.18% and 99.36%, respectively. The accuracy of the HCB was 94.60%, and the extraction effect of the CW was not perfect.
The accuracy of the diameter at different positions in the stem analysis was between 85% and 98%. The extraction effect at RH10 was the best (98.43%), but it was the worst at RH90 (85.51%). The variation trend of the diameter extraction accuracy is as follows: from RH2 to RH10, it increases with height and reaches a maximum at RH10; then, it gradually decreases as the height increases. In this study, the accuracy at RH80 was higher than that at RH75, which may be due to the shaking of the tree crown, and the multi-station point cloud data could not be completely matched (as shown in Figure 2), which caused a large visual error. However, at RH80, the stem of each station was separated from the others, and the visual error was relatively reduced. At RH90, the point cloud data were missing, the density was reduced, and the error was increased.   Table 4 shows the fit results of the taper equation. Both datasets have the highest fitting accuracy of model (6). The fitting effect of the TLS data in models (1)-(3) was better than the measured data; however, models (4)-(6) yielded the opposite results. In general, the difference between the TLS data and stem analysis data was very small; i.e., R 2 > 0.91. The difference between the optimal models was only 0.0024, and the Bias and RMSE were also very close.  Table 4 shows the fit results of the taper equation. Both datasets have the highest fitting accuracy of model (6). The fitting effect of the TLS data in models (1)-(3) was better than the measured data; however, models (4)-(6) yielded the opposite results. In general, the difference between the TLS data and stem analysis data was very small; i.e., R 2 > 0.91. The difference between the optimal models was only 0.0024, and the Bias and RMSE were also very close. To more comprehensively evaluate the fitting effect of the TLS data, Figures 3 and 4 plot the residual graphs and trend lines of the six models of both datasets, respectively. The figure shows that both datasets have the high equal variance and unbiased characteristic of model (6), and the residual distribution is more uniform. To more comprehensively evaluate the fitting effect of the TLS data, Figures 3 and 4 plot the residual graphs and trend lines of the six models of both datasets, respectively. The figure shows that both datasets have the high equal variance and unbiased characteristic of model (6), and the residual distribution is more uniform.       (5), the R 2 of the TLS data of the other models as greater than the measured data. The prediction accuracy of the TLS data was more than 97%, and the testing results of the six models were consistent with the fitting effect. Both datasets had the largest R 2 in model (6), the smallest    (5), the R 2 of the TLS data of the other models as greater than the measured data. The prediction accuracy of the TLS data was more than 97%, and the testing results of the six models were consistent with the fitting effect. Both datasets had the largest R 2 in model (6), the smallest Bias and RMSE, and the highest accuracy. Therefore, the model (6) among these was selected and used for the next analysis. To more clearly describe the effect of the taper equation established by the TLS data, the average values of the DBH and TH extracted from the sample plot were used as the DBH (17.20 cm) and TH (18.64 m) of the simulated tree, and model (6) was used to draw the stem curve of both datasets. As shown in Figure 5, the taper model of TLS basically coincides with the stem analysis. According to the fitting and testing results of stem analysis, the effects of the TLS data are consistent with the stem analysis data, with only slight differences in statistical indicators. Therefore, we believe that TLS data can be used instead of the data from felled trees to establish the taper equation.

Fitting and Testing of the Taper Models with Different Relative Heights
Since the scanning effect of TLS is easily blocked by branches and leaves, it is difficult to obtain detailed structural information on the upper part of the tree. When the diameter is extracted to RH50, the accuracy is 96.78% and decreases with height. The accuracy of the TLS data will affect the fitting effect of the taper model. Therefore, according to the extraction accuracy mentioned above, this study used diameters from 0 m to relative heights of RH50, RH60, RH70, RH75, RH80, and RH90 to establish the taper equation. The optimal basic taper equation (model 6: Kozak (2002)-II) mentioned above was adopted, and the same method was used for fitting and testing. The results are shown in Table 6.

Fitting and Testing of the Taper Models with Different Relative Heights
Since the scanning effect of TLS is easily blocked by branches and leaves, it is difficult to obtain detailed structural information on the upper part of the tree. When the diameter is extracted to RH50, the accuracy is 96.78% and decreases with height. The accuracy of the TLS data will affect the fitting effect of the taper model. Therefore, according to the extraction accuracy mentioned above, this study used diameters from 0 m to relative heights of RH50, RH60, RH70, RH75, RH80, and RH90 to establish the taper equation. The optimal basic taper equation (model 6: Kozak (2002)-II) mentioned above was adopted, and the same method was used for fitting and testing. The results are shown in Table 6. The taper model of 0 m-RH60 has the highest accuracy (99.16%), which is 0.15% higher than that obtained using all the data. At 0 m-RH50, the accuracy is high, but there are fewer data representing the shape of the upper stem. At 0 m-RH70, RH75, RH80, and RH90, the accuracy of the model is consistent with the diameter extraction, and both gradually decrease with height.
Based on the accuracy of diameter extraction and model fitting, RH70 is considered to be the best height for TLS to obtain the stem diameter. The extraction accuracy of RH70 is slightly different from RH60 but is much higher than that of RH75 and RH80. As more data are collected from RH70 to RH90, the lower the accuracy of the model is (99.12% at RH70). Therefore, RH75, RH80, RH90 were discarded and the data having 0 m-RH70 diameters were selected for further analyses.  Table 7. The t-test shows that the model parameter estimates are all significant (p < 0.0001), and the fitting statistics (R 2 = 0.9665) indicate that the model could well describe the variation trend of the stem shape of larch.  Table 8 lists some of the models with a good fitting of random parameters. It can be seen from the table that the mixed model with random effects is superior to the basic model. When a random parameter a 8 is introduced, the fitting effect of the mixed model is better; a 6 and a 9 , as the combination of two random parameters, have a better fitting effect. The likelihood ratio test shows that the models are significantly different from each other. A comprehensive analysis considering the tree effect with two random parameter combinations (a 6 and a 9 ) is the best mixed model.

Heteroscedasticity and Autocorrelation of Mixed Model Errors
In this study, the heteroscedasticity of the error was judged by the residual distribution diagram. The residual in Figure 6 shows the trend of the random distribution, which indicates that the errors of the mixed model have no heteroscedasticity. First-order autoregressive structure AR (1) is used to solve the autocorrelation problem of diameters measured on the same tree. Table 9 shows that when AR (1) is added, the three evaluation indexes and likelihood ratio confirmed that the AR (1) significantly improved the fitting accuracy of the mixed-effect model.   Table 10 lists the parameter values and fitting statistics of the basic model and the mixed model with the tree effect. The results show that the R 2 of the mixed model is greater than that of the basic model, and the Bias and RMSE are both smaller, which indicates that the fitting effect of the mixed model is better than that of the basic model.    Table 10 lists the parameter values and fitting statistics of the basic model and the mixed model with the tree effect. The results show that the R 2 of the mixed model is greater than that of the basic model, and the Bias and RMSE are both smaller, which indicates that the fitting effect of the mixed model is better than that of the basic model. Note: σ 2 i and σ 2 ij are the variance and covariance parameters for the tree effect; σ is the residual variance; and ρ is an autoregressive parameter, ρ > 0.

Mixed Model Test
The prediction accuracy of the random effect mixed model reached 99.72%, which is higher than that of the basic model and fixed effect model, and the bias was reduced by 0.0249 cm (Table 11). This result indicates that the taper model based on the TLS data can be widely used in forestry. The accuracy (Table 12) of the taper equation in deriving the diameter can reach more than 85% (except RH90). Consistent with the variation trend of the TLS extraction accuracy, the lower accuracy of the upper stem is only 74.22% at RH90. In forest inventory assessments, tree crowns and treetops are generally not used for timber production. Therefore, the data obtained using TLS in this study have important practical value.

Volume Estimation by the Mixed Model
The volume predicted by the two methods is very close to the reference value. Method 2 provides better performance (Table 13), i.e., the method of nondestructive TLS measurements is suitable for deriving the taper function and stem volume.

Discussion
The conventional construction of the taper equation is performed by obtaining data by cutting down trees, which is a waste of time and energy. This research provides a new nondestructive TLS scanning method, which not only provides Woodsfield data for verification but also provides large-scale plot data for modeling. Compared with Sun's study (applying 198 trees as samples in 8 plots, 16 of which were subjected for stem analysis) [32], this study involved more data (a total of 580 trees in 9 plots, and another 30 were subjected to stem analysis), and this study therefore better represents the actual state. Most studies on the diameter of TLS acquisition are limited to stems with DBHs, while there are few studies on upper TLS acquisition. Xi et al. [21] reported that the occlusion of the crown introduces biases above 12 m in stem form extraction. In this study, the diameters of the relative positions were extracted, which not only reduced the influence of the tree height but also derived the diameters through the taper model.
According to the diameter extraction of different relative heights in the stem analysis, the accuracy of the entire tree decreases with height. Notably, the accuracy at RH75 is less than that at RH80. The study showed that during TLS scanning, due to the influence of wind blowing on the tree or the error in the registration of the station, there are angles or ghosts in the upper part of the stem at different stations. There will be errors in the registration of different stations of point cloud data, which is shown in Figure 2a, and tree heights in this study ranged from RH75 to RH80. The visual error, shown in blue in Figure 2b, may overestimate or underestimate the diameter. At RH75, the tree point clouds are mixed together, so it is not easy to interpret the stem. At RH80, the stem point clouds are separated, which is easy to interpret and has little visual error. At 0.9H, the point cloud density decreases, the stem diameter is difficult to interpret, and the extraction error increases. Beam divergence of the scanner is the key determinate of the increasing error with height. To overcome this problem, we used multiple scans of trees from different perspectives. Xi et al. [21] scanned plots using Leica HDS6100 by the Finnish Geodetic Institute (FGI), and the results showed that the DBH extraction accuracy was 0.97 (r 2 ) and 0.90 cm (RMSE). In our study, the DBH extraction accuracy was 0.9958 (R 2 ) and 0.3726 cm (RMSE). Li Dan et al. [33] used Riegl VZ-400 to scan 5 plots, and the average RMSE of DBH and TH extracted automatically were 2.53 cm and 3.80 m, respectively (lower than ours). In addition, the stand density, crown density, tree height, branch growth condition, understory vegetation, and number and distance of stations will all affect the accuracy.
Internationally, the tree detection rate, stem position, DBH, TH, wood volume, biomass and other forest measurement processes of TLS have reached a high degree of automation [34]. Compared with traditional methods, the DBH and TH are automatically obtained by algorithms, such as Cabo C, where 85% of the DBH deviations are lower than 1 cm and 92% of the TH differences are less than 0.5 m [35]. This study used manual extraction, except for 0 m, and all diameter errors were less than 1 cm. With regard to the quality of the data, the accuracy in this study was higher, but the processing time was longer and the workload was greater. Regardless, the lower portions of the trees are obscured by understory vegetation, as well as other trees and branches, or the lack of scanning points on the treetops.
According to fitting analysis of the three types of taper equations, it was concluded that the variable-exponential taper equation has the highest accuracy, and the segmented taper equation also has a higher accuracy, which is consistent with previous studies. Fitting the basic taper equation, the results indicate that RH70 is the maximum height for TLS to obtain diameter data. In our view, this position has higher accuracy and can reduce the processing time of redundant data.
Relevant studies have shown the predictive ability of the random effect mixed model preference at tree levels [36]; therefore, this study only considered the tree effect but not the plot effect. A stem taper is susceptible to the tree crown characteristics. Whether it is a segmented taper equation [12] or a variable-exponential taper equation, the mixed model with the tree effect can better manage the hierarchical structure of the data. Another innovation of this study is that TLS was entirely used to obtain the diameters, and TH was used to establish a taper equation mixed model, which was obtained without cutting down trees. The RD1000 was also adopted for non-destructive tree measurement, as it had the largest mutation when the relative height was 0.64-0.8, similar to the results of our previous study (RH75) [37]. Sun et al. [32] used a modified Schumacher equation to fit the Populus L. taper equation of the TLS, which yielded a fitting R 2 = 0.96, which is similar to that of the basic model in our study but lower than that of the mixed model (R 2 = 0.98). Other studies have recommended that, as a variable for calculating volume, the cross-section area exhibits better performance than the diameter [38]. The prediction quality of the taper and volume model can also be calibrated with an additional diameter measured between 40% and 60% of the total tree height [36,39,40]. When the diameter of the calibration measurement is 7 m, the average deviation of the merchantable volume is only 0.63% [41].

Conclusions
In this study, TLS was used to obtain forest parameters, a small amount of field data was used to verify its accuracy, and a mixed model was established to improve the prediction effect of the taper equation. The results show that nondestructive TLS measurements are suitable for deriving stem diameters and trunk volumes; they are feasible for forest investigations.
The optimal height of diameter extraction is RH70; this accuracy meets the requirements of forestry surveys, and also reduces much of the point cloud data processing work. If the extracted height is too high, the diameter accuracy of the upper tree will be low, which will cause a significant amount of redundant data. Conversely, if the extracted height is too low, the amount of data acquired will be too low, which will affect the stem curve. In addition, the amount and accuracy of the data will affect the fitting of the taper model, which in turn affects the simulation of the stem curve and volume estimation.
The taper equation model based on both datasets has a better fitting effect and higher accuracy (model (6): R 2 > 0.97, P > 98.94%), so the optimal basic taper equation was selected to establish the mixed model with the tree effect. The accuracy of the mixed model based on the tree effect reached 99.72%. The mixed model can accurately predict the diameter of the upper stem, which resolved the problem of the occlusion by the crown. The diameter accuracy derived from the mixed model reached 99%, and the bias was <1.5 cm; it also had high accuracy in estimating the volume of the stem.
The taper equation established by TLS data has typical representativeness and good practicability. In future field measuring work, this method can be used to conduct highaccuracy forest surveys without cutting trees.