Accounting for terrain in neighborhood competition: fused UAV and backpack LiDAR for fine-scale volume estimation in Chinese fir (Cunninghamia lanceolata) plantations

1 Introduction

Chinese fir (Cunninghamia lanceolata) is an important fast-growing timber species in southern China. Its plantations span 10.9 million hectares and possess a growing stock volume of 993.37 million cubic metres, representing approximately 19% of the national plantation area and 26% of the national plantation growing stock volume (Jiang et al., 2022). As individual-tree volume is the basic unit used to calculate stand growing stock volume, accurate estimation of individual-tree and stand growing stock volume under operational forest management conditions is important for assessing silvicultural effectiveness, optimizing stand density and spatial structure, and supporting fine-scale forest carbon accounting (Omoniyi and Sims, 2024; Yang Z. et al., 2025).

Although traditional field-based growing stock volume surveys remain a vital method for forest resource inventory, they are often costly and labor-intensive, with their efficiency significantly hampered in complex terrain and closed-canopy stands, where structural information may be incompletely captured (Omoniyi and Sims, 2024). To address these limitations in the acquisition of three-dimensional structural data, light detection and ranging (LiDAR) technology has increasingly been utilized for the retrieval of forest structural parameters and the estimation of forest attributes (Coops et al., 2021; Maeda et al., 2025). Unmanned aerial vehicle light detection and ranging (UAV LiDAR) can rapidly cover relatively large areas and effectively retrieve upper-canopy attributes, but its ability to capture understory and stem information is often limited by canopy occlusion, particularly in stands with dense foliage or high canopy closure (Fraser et al., 2025; Li Y. et al., 2025). In contrast, backpack laser scanning (BLS), as a form of mobile laser scanning (MLS), enables continuous field data collection and provides detailed stem information more rapidly than static terrestrial laser scanning (TLS) surveys (Xu Z. et al., 2023; Li Y. et al., 2025). Integrating UAV LiDAR and backpack LiDAR data within a common coordinate system and terrain framework enables the acquisition of more comprehensive structural information from both the upper canopy and the understory. This integration has been shown to improve the accuracy and stability of forest attribute estimation, including growing stock volume (Fekry et al., 2022; Balestra et al., 2024; Ali et al., 2025). However, despite the enhanced structural information provided by multi-source LiDAR, many studies continue to oversimplify neighborhood relationships when modeling tree competition (Huang et al., 2024; Rieder et al., 2024).

In numerous studies using LiDAR or multi-source LiDAR data to estimate forest attributes and examine tree competition, stand structure is often simplified through two-dimensional (2D) or two-and-a-half-dimensional (2.5D) analyses. For instance, canopy height models (CHMs) are frequently employed for tree height extraction and volume estimation (Daryaei et al., 2025). While these methods are valuable for operational forest inventories, they may not fully capture the within-stand spatial arrangement and vertical differentiation. This limitation is particularly evident in the quantification of neighborhood competition, as tree competition is commonly quantified using horizontal distance alone, disregarding elevation differences and terrain conditions (Rieder et al., 2024; Ronoud et al., 2024). Consequently, variables associated with competition may inadequately reflect the relationships inherent in neighborhood competition within stands exhibiting significant crown overlap and marked topographic variation. However, whether three-dimensional competition metrics consistently improve tree-level prediction remains uncertain, because their contribution may depend on competitor representation, focal-tree attributes, data quality, and modeling strategy (Rieder et al., 2024; Ronoud et al., 2024). Therefore, the contribution of competition indices of varying dimensions to the accuracy of individual-tree volume estimation still necessitates further investigation.

Recent point-cloud-based research has tackled the above limitations by investigating 3D forest structure and tree competition across various forest types. Pont et al. (2021) found that the integration of inter-tree competition metrics derived from airborne laser scanning into spatial models of radiata pine accounted for previously unexplained variation in tree traits, including stem volume, and reduced residual variance. Irwin et al. (2025) further found that, in mixed stands of Douglas-fir, interior spruce, and lodgepole pine, drone- LiDAR-derived competition metrics were significant predictors of basal area growth when assessed individually, whereas focal-tree crown volume provided the strongest predictive information. Complementing these point-cloud-based studies, Long et al. (2020) demonstrated that, for ring-cupped oak trees, accounting for slope and aspect and incorporating an improved Hegyi competition index significantly improved individual-tree height prediction, underscoring the importance of terrain information in competition-related modeling. However, in studies of Chinese fir plantations, empirical research remains limited on the construction of digital elevation model (DEM)-constrained three-dimensional competition indices based on integrated UAV LiDAR and BLS data and on the evaluation of their contribution to individual-tree volume prediction across different stand density gradients. A small body of literature has compared different formulations of inter-tree competition metrics. For example, Ronoud et al. (2022) contrasted stem-derived and crown-based competition indices derived from terrestrial laser scanning data in Scots pine stands. Nevertheless, few controlled comparative analyses have held the model structure and focal-tree traits constant while varying only the dimensionality of the competition indices.

Accordingly, the overall objective of this study was to determine whether DEM-constrained three-dimensional competition indices could improve the accuracy of individual-tree volume estimation, compared with the conventional two-dimensional competition indices, in Chinese fir plantations. To the best of our knowledge, this study is the first to integrate the complementary canopy and stem information provided by UAV LiDAR and backpack LiDAR with DEM-derived true three-dimensional inter-tree distances and slope-weighted competition indices and to apply this integrated framework to individual-tree volume prediction in Chinese fir plantations. The specific objectives were to: (1) develop a terrain-aware framework for deriving DEM-constrained three-dimensional competition indices from fused UAV LiDAR and BLS data; (2) determine whether three-dimensional competition indices improve the accuracy of individual-tree volume prediction compared with conventional two-dimensional competition indices, and whether their effects vary across stand-density levels, input organizations, and modeling approaches; and (3) use Shapley additive explanations (SHAP) to quantify the relative contributions of focal-tree attributes, neighbor attributes, and competition variables to individual-tree volume prediction.

2 Materials and methods

2.1 Study area

The study area is located in the Fengshushan Forest Farm of Jingdezhen City, Jiangxi Province, China (28°44′–29°56′N, 116°57′–117°42′E; Figure 1), which belongs to the subtropical evergreen forest distribution areas in the southern China. Among all planted evergreen coniferous tree species, Chinese fir accounts for the largest plantation area, approximately 31.0% of the total farm area (Zeng et al., 2015). In Chinese fir plantations, the shrub layer is mainly composed of Chinese fringe flower (Loropetalum chinense) and Japanese maesa (Maesa japonica), and the herb layer mainly includes Japanese chain fern (Woodwardia japonica) and Old World forked fern (Dicranopteris pedata). The Fengshushan Forest Farm is located in a transitional landscape between the foothills of the Huangshan and Huaiyu Mountains and the Poyang Lake Plain, exhibiting significant topographic relief. The 12 sample plots were situated at elevations ranging from 121.9 to 180.4 m above sea level. The soils in this region are predominantly red soils derived from granite parent material. The area experiences a humid mid-subtropical monsoon climate characterized by warm and moist conditions, with a mean annual air temperature of 17.1 °C and an average annual precipitation of approximately 1,805 mm (Yin et al., 2024). Precipitation in the study area is strongly seasonal, with December falling within the relatively low-precipitation period of the year. These relatively dry conditions were favorable for field surveys and facilitated the acquisition of high-quality LiDAR point-cloud data (Kang et al., 2020).

2.2 Field inventory and single-tree volume

The Chinese fir experiment was established in the spring of 1999 with three planting-density treatments. Each planting density was replicated five times in a randomized block design. Each plot comprised a rectangular area of 30 m × 20 m (600 m²) surrounded by a one-row buffer zone consisting of similarly treated trees. Of the 15 original plots, two plots from the original low-density treatment and one plot from the original medium-density treatment were excluded because they had been converted from pure Chinese fir stands into Chinese fir–broadleaved mixed forests, leaving 12 25-year-old pure Chinese fir plots. Thus, field inventory data were collected in December 2024 from all 1,453 trees in the 12 retained plots. Based on the tree counts recorded in this inventory, stand density was calculated and expressed as stems ha⁻¹. According to the thresholds of Liu and Tong (1980), the retained plots were classified as low-, medium-, and high-density plots (< 1,667, 1,667–2,000, and > 2,000 stems ha⁻¹, respectively), yielding three, four, and five plots in the respective classes. A GNSS real-time kinematic (RTK) receiver (South P30, South GNSS Navigation Instrument Co., Ltd., China) was used to record the coordinates of four corners of each plot and stem coordinates of all individual trees within the plots. Individual-tree measurements included diameter at breast height (DBH, cm), tree height (H, m), and crown width (CW, m). DBH was measured at 1.3 m above ground. Tree height was measured using a height pole, and crown width was measured with a tape along two orthogonal directions and calculated as the mean of the two measurements. Plot-level summary statistics are provided in Table 1. In this study, individual-tree volume, calculated from field-measured diameter at breast height (DBH) and tree height (H) using an empirical two-variable volume equation for Chinese fir, was used as the reference response variable (Liu and Tong, 1980).where V o b s e r v a t i o n is the equation-derived reference individual-tree volume (m³), D B Hdenotes diameter at breast height (cm), and Hdenotes tree height (m).

V o b s e r v a t i o n = 0.000058777042 D B H 1.9699831 · H 0.89646157 (1)

Table 1. Summary of stand characteristics measured during the December 2024 field inventory across 12 Chinese fir plantation plots.

2.3 LiDAR acquisition and point cloud processing

LiDAR data were acquired in December 2024 using a DJI Matrice 300 RTK equipped with a Zenmuse L1 sensor for unmanned aerial vehicle light detection and ranging acquisition and a GreenValley DGC50H system for backpack laser scanning (Xu D. et al., 2023; Xu Z. et al., 2023; Gao et al., 2025). UAV LiDAR data were acquired using an automatically planned single-grid flight route for LiDAR mapping in the DJI Pilot platform (Ye et al., 2025; Jarahizadeh and Salehi, 2026). No ground control points were deployed during aerial surveys, and RTK positioning corrections were obtained via 5G network connection to a continuously operating reference station (CORS) service (Gómez-Gutiérrez et al., 2024). During UAV LiDAR data acquisition, two-dimensional visible-light images were simultaneously captured by the integrated RGB camera of the Zenmuse L1 sensor; however, these images were not processed or utilized in the subsequent analyses of this study. For UAV LiDAR acquisition, a terrain-following mode was utilized, with a flight altitude of 100 m, a forward overlap of 70%, a side overlap of 80%, and an average point density of approximately 260 points m⁻². For backpack laser scanning acquisition, each plot was scanned along an S-shaped trajectory at an approximate walking speed of 0.5 m s⁻¹ to enhance the capture of stems and understory structure. The average point density of the backpack LiDAR data across the 12 plots was approximately 27,096 points m⁻².

The denoised UAV LiDAR and BLS point clouds were first transformed to the China Geodetic Coordinate System 2000 (CGCS2000) using a 3° Gauss–Krüger projection with a central meridian of 117°E (Yang, 2009). To eliminate pose deviations between UAV LiDAR and BLS point clouds, overlapping mid-canopy point slices within a height range of 2–10 m above ground were extracted for preliminary coarse alignment, followed by fine registration using the iterative closest point (ICP) algorithm (Xu et al., 2021; Shao et al., 2022; Zhang et al., 2025). For ICP refinement, point clouds were downsampled to a voxel size of 0.10 m, and a correspondence search radius of 0.30 m was employed, with 2% of correspondences trimmed as outliers (Shao et al., 2022; Zhang et al., 2025). Within each 0.10 m voxel, redundant points were represented by a single average voxel point (Liu H. et al., 2020; Lyu et al., 2024). The 0.30 m search radius was selected as three times the voxel size, which is consistent with recommendations that the matching distance should be approximately two to four times the sampling voxel size. In instances where automated alignment did not reach stable residual patterns, only the initial transformation was manually adjusted based on visually identifiable correspondences of stems and canopies, after which ICP refinement was repeated using the same parameters. After accurate registration, the UAV LiDAR and BLS point clouds were fused into a single integrated point cloud for each plot, which was used for subsequent preprocessing, individual-tree segmentation, and forest attribute extraction (Xu D.et al., 2023; Balestra et al., 2024). After denoising, the average point density of the fused UAV LiDAR and BLS point clouds was approximately 27,343 points m⁻².

The fused point cloud was denoised and height-normalized using the ground surface (Li Y. et al., 2025). Stem points were extracted using Density-Based Spatial Clustering of Applications with Noise (DBSCAN) clustering in conjunction with the Circular Slice Projection (CSP) (Zhang et al., 2022). DBH was measured in LiDAR360 by fitting a circle to stem points within a cylindrical region extending from 1.25 to 1.35 m above ground, with a horizontal radius of 0.15 m. The fitted diameter was recorded as DBH, and the circle center was used as the seed point for individual-tree segmentation (Xu Z. et al., 2023; Li Y. et al., 2025). All seed points were visually inspected and corrected as necessary. Individual-tree segmentation was performed using the Point Cloud Segmentation (PCS) method, with stem seeds serving as inputs, implemented in LiDAR360 (GreenValley International) (Deng et al., 2024; Li Y. et al., 2025). For each segmented tree, H was obtained from the maximum height-normalized point, and CW was calculated as the mean of the crown diameters along two orthogonal horizontal directions. Segmented trees were matched to field trees based on the spatial correspondence between fitted stem-center coordinates and RTK-measured stem positions. Segmentation performance was evaluated using recall, precision, and F1-score at the detection stage. The accuracy of LiDAR-derived DBH, H, and CW was evaluated using linear regression metrics (R², RMSE, and MAE) by comparison with the corresponding field measurements for each tree (Li Y. et al., 2025; Zhang et al., 2025).

2.4 DEM-constrained neighborhoods, distances, and competition indices

Based on planar coordinates (x, y), and given its adaptive neighbor delineation mechanism independent of predefined threshold settings and its robustness across varying stand structures, Voronoi tessellation and Delaunay triangulation were employed to identify neighbors (Figure 2) (Li et al., 2021; Zhou et al., 2022). Voronoi tessellation partitions a planar domain into a set of polygonal cells based on the nearest-neighbor rule, where any point within a specific cell is geographically closer to its corresponding seed point than to any other seed point. As the geometric dual of Voronoi tessellation, Delaunay triangulation connects seed points whose Voronoi cells share an edge and satisfies the empty-circle criterion, whereby no additional seed point lies within the circumcircle of any Delaunay triangle (Voronoi, 1908; Delaunay, 1934). For each focal tree, all trees connected to it in the Delaunay triangulation were treated as neighboring trees. Delaunay triangulation restricted neighbors to locally adjacent trees, while inverse-distance weighting reduced the influence of more distant neighbors. A 0.1 m DEM was generated from the ground-classified points of the fused point cloud using triangular irregular network (TIN) interpolation in LiDAR360. The grid size was selected based on the density and spatial continuity of the ground-classified points within the plots and the need to retain fine-scale terrain differences between neighboring trees in plots with marked topographic variation. The DEM was then used to obtain three-dimensional coordinates (x, y, z) and calculate the three-dimensional distance and slope angle between focal trees and their neighbors (Figure 3) (Rieder et al., 2024). To mitigate boundary truncation effects during Voronoi–Delaunay neighborhood construction, a 2 m buffer was established outside each plot, following previous forest spatial-structure studies that adopted the same buffer width; trees within the buffer were included only as potential neighboring trees of focal trees within the original plot boundary and were not treated as focal trees (Pommerening and Stoyan, 2006; Liu et al., 2023; Wang et al., 2023).

The two-dimensional and three-dimensional inter-tree distances and the slope angle were calculated as follows (Equations 2–4):where d 2 D , i jand d 3 D , i jare the two-dimensional and three-dimensional inter-tree distances (m), respectively, and θ i jis the terrain slope angle (radians) along the horizontal line connecting focal tree i and neighboring tree j.

d 2 D , i j = ( x j − x i ) 2 + ( y j − y i ) 2 (2)

d 3 D , i j = ( x j − x i ) 2 + ( y j − y i ) 2 + ( z j − z i ) 2 (3)

θ i j = arctan ( | z j − z i | ( x j − x i ) 2 + ( y j − y i ) 2 ) (4)

Two input organizations were considered (Figure 4). (1) The competition-pair retained neighbor DBH, H, and CW, inter-tree distances, slope angles, and the single-neighbor contribution terms within the Hegyi, Arney, and CCI summations for each focal–neighbor pair (1v1). (2) The structural-unit averaged these variables across all neighbors of each focal tree to produce one sample per focal tree (1vN). A total of six distance-based competition indices, three conventional two-dimensional indices and their three DEM-constrained three-dimensional counterparts, were constructed based on the Hegyi index (DBH-distance), the Arney index (H-distance), and the Crown Competition Index (CCI; CW-distance) (Corona and Ferrara, 1989; Wang et al., 2024). To avoid zero denominators, a constant ε = 10 − 3m was added to all denominators (Firmino et al., 2023). For each of the three index types, both a conventional two-dimensional version and a DEM-constrained three-dimensional version were calculated. Within a given index and dimensional formulation, larger values indicate greater estimated competitive pressure exerted by neighboring trees on the focal tree.

Two-dimensional indices used d 2 D , i j (Equations 5–7):

H e g y i i D B H − 2 D = ∑ j ∈ N ( i ) D B H j D B H i d 2 D , i j + ε (5)

A r n e y i H − 2 D = ∑ j ∈ N ( i ) H j H i d 2 D , i j + ε (6)

C C I i C W − 2 D = ∑ j ∈ N ( i ) C W j C W i d 2 D , i j + ε (7)

Three-dimensional indices replaced d 2 D , i j with d 3 D , i j and incorporated DEM-derived slope effects through a slope-weighting term (Equations 8–11):

ω θ = exp ( − k θ θ i j ) , k θ > 0 (8)

H e g y i i D B H − 3 D = ∑ j ∈ N ( i ) D B H j D B H i d 3 D , i j + ε · ω θ (9)

A r n e y i H − 3 D = ∑ j ∈ N ( i ) H j H i d 3 D , i j + ε · ω θ (10)

C C I i C W − 3 D = ∑ j ∈ N ( i ) C W j C W i d 3 D , i j + ε · ω θ (11)

where k θwas set to a fixed constant value (k θ = 0.02 rad − 1) in our implementation. This avoids fold-specific data-driven tuning and keeps the evaluation focused on isolating the dimensionality effect of the competition indices under consistent model settings. The slope weight ω θ was computed from the geometry prior to model training and then used consistently across all cross-validation folds.

2.5 Input variables and modeling framework

This study adopted five machine learning models and three deep learning models, including Random Forest (RF) (Breiman, 2001), Extra Trees (ET) (Geurts et al., 2006), Extreme Gradient Boosting (XGBoost) (Chen and Guestrin, 2016), Light Gradient Boosting Machine (LightGBM) (Ke et al., 2017), histogram-based gradient boosting (HistGBR) (Guryanov, 2019), multilayer perceptron (MLP) (Murtagh, 1991), Transformer (Vaswani et al., 2017), and graph convolutional network (GCN) (Kipf and Welling, 2016). These algorithms are representative and widely used data-driven approaches for remote sensing-based forest attribute retrieval. The set of candidate models was deliberately selected to encompass mainstream algorithm families with robust and validated predictive performance in forest biometric modeling. Overall, these eight models have been extensively utilized to estimate core forest productivity indicators, such as stand volume, biomass, and carbon storage, in previous forestry studies (Seely et al., 2023; Omoniyi and Sims, 2024; Liu et al., 2025; Miao et al., 2025; Yang L. et al., 2025; Qian et al., 2026). Input variables included LiDAR-derived focal-tree growth attributes (DBH, H, and CW), LiDAR-derived neighbor attributes (neighbor DBH, H, and CW), neighbor count, distance metrics, slope angle, and competition indices. These variables were used to predict individual-tree volume (V p r e d), with validation accuracy assessed by comparing V p r e d to the reference individual-tree volume (V o b s e r v a t i o n) calculated from field-measured DBH and H. Input variables were organized according to the competition-pair and structural-unit configurations described above.

All numeric predictors were standardized using z-score normalization, including focal-tree DBH, H, and CW; neighbor DBH, H, and CW or their mean values; distance and angle variables included in the corresponding configuration; Hegyi, Arney, and CCI competition indices; and neighbor count in the structural-unit configuration. To prevent information leakage, the parameters for standardization were computed from the training set and subsequently applied to the validation split within each fold (Little et al., 2017; Xu and Goodacre, 2018; Wilimitis and Walsh, 2023; Apicella et al., 2025). The performance of the model was evaluated using five-fold cross-validation grouped by focal tree, ensuring that all samples associated with the same focal tree were assigned to the same fold (Little et al., 2017). In each fold, approximately 80% of focal trees were used for training and the remaining 20% for validation (Xu and Goodacre, 2018). Results were summarized using the mean ± standard deviation (SD) of the coefficient of determination (R²), root mean square error (RMSE), and mean absolute error (MAE) across the five validation folds (Wilimitis and Walsh, 2023). For feature attribution, SHAP values were calculated for the held-out validation samples in each cross-validation fold using TreeExplainer implemented in the SHAP Python package (Lundberg et al., 2020; Scheda and Diciotti, 2022). For each predictor, the mean absolute SHAP value was calculated within each fold, averaged across the five folds, and normalized across all predictors to obtain its relative contribution (%) (Dieste et al., 2024; Schneider et al., 2026). Positive and negative SHAP values indicate increases and decreases in the predicted individual-tree volume relative to the model baseline, respectively, whereas absolute SHAP values indicate the magnitude of the predictor contribution (Lundberg and Lee, 2017). For the reported experiments, each model was executed with fixed hyperparameter settings within each configuration to maintain a controlled comparison across input structures and competition dimensionality. All models were implemented in a Python 3.8 environment.

3 Results

3.1 Tree growth parameter extraction

After excluding the largest 2% of correspondence residuals as potential mismatches, the plot-level registration RMSE was approximately 0.06 m for most plots, with the highest value of 0.072 m observed in plot H-5. These residuals were sufficiently small to support subsequent seed-based individual-tree segmentation and the calculation of inter-tree distances, slope angles, neighbor counts, and slope-weighted competition indices.

Individual-tree segmentation based on the fused point cloud achieved a mean recall of 0.91 ± 0.01, a mean precision of 0.93 ± 0.01, and a mean F1-score of 0.92 ± 0.01, where the values are reported as the mean ± standard deviation across the 12 plots. As shown in Table 2, LiDAR-derived DBH, H, and CW showed good agreement with field measurements. When pooled across all trees from all plots, the fitted relationships were DBH: y = 0.85 x + 2.52(R² = 0.91, RMSE = 1.52 cm, MAE = 1.16 cm), H: y = 0.95 x + 0.81 (R² = 0.89, RMSE = 1.28 m, MAE = 1.00 m), and CW: y = 0.84 x − 0.68 (R² = 0.74, RMSE = 0.45 m, MAE = 0.37 m), where (y) represents the field-measured value from the sample plots, and (x) represents the corresponding value extracted from the fused point cloud for DBH, H, or CW. Based on the pooled R² values, DBH showed the highest agreement with field measurements, followed closely by H, whereas CW showed lower agreement.

Table 2. Accuracy assessment of LiDAR-derived structural parameters (DBH, H, and CW) compared with field measurements across 12 Chinese fir plots.

3.2 Distribution of the competition indices

After incorporating DEM-based slope correction, the DEM-constrained three-dimensional competition indices had smaller values than the corresponding two-dimensional competition indices across both input organizations (Figures 5, 6). In the competition-pair configuration (1v1; Figure 5), the overall medians of the three-dimensional Hegyi, Arney, and CCI indices were 0.466, 0.353, and 0.406, respectively, representing reductions of 12.5, 11.6, and 11.6% relative to the corresponding two-dimensional indices. Similarly, in the structural-unit configuration (1vN; Figure 6), the overall medians of the three-dimensional Hegyi, Arney, and CCI indices were also lower than their two-dimensional indices, with decreases of 11.8, 11.9, and 11.8%, respectively. Across stand-density classes, high-density stands generally showed higher median index values, while the density-class ordering varied among indices and input configurations. Within each stand-density class, input organization, and dimensional formulation, the Hegyi index consistently displayed the highest median value among the three competition indices.

3.3 Model performance comparison

Across input organizations and stand-density classes, RF showed the most consistent improvement from 2D to 3D inputs, with the largest gains under the 3D competition-pair configuration (Tables 3–6). Compared with the 2D competition-pair configuration, the 3D competition-pair configuration increased R² from 0.73 to 0.87 in low-density stands, from 0.56 to 0.76 in medium-density stands, and from 0.81 to 0.90 in high-density stands. Simultaneously, RMSE decreased from 0.0253, 0.0441, and 0.0465 to 0.0182, 0.0264, and 0.0266 m³ tree⁻¹ in low-, medium-, and high-density stands, respectively. At the overall scale, R² increased from 0.74 to 0.85 and RMSE decreased from 0.0434 to 0.0346 m³ tree⁻¹, indicating that RF achieved higher individual-tree volume prediction accuracy with the 3D competition-pair configuration than with the corresponding 2D configuration.

Table 3. Model performance (mean ± SD of R², RMSE, and MAE) across input settings (2D competition-pair, 2D structural-unit, 3D competition-pair, and 3D structural-unit) and stand-density classes (Low, Medium, High, and Overall): RF and LightGBM.

Table 4. Model performance (mean ± SD of R², RMSE, and MAE) across input settings (2D competition-pair, 2D structural-unit, 3D competition-pair, and 3D structural-unit) and stand-density classes (Low, Medium, High, and Overall): Extra Trees and XGBoost.

Table 5. Model performance (mean ± SD of R², RMSE, and MAE) across input settings (2D competition-pair, 2D structural-unit, 3D competition-pair, and 3D structural-unit) and stand-density classes (low, medium, high, and overall): HistGBR and MLP.

Table 6. Model performance (mean ± SD of R², RMSE, and MAE) across input settings (2D competition-pair, 2D structural-unit, 3D competition-pair, and 3D structural-unit) and stand-density classes (low, medium, high, and overall): transformer and GCN.

In contrast, other tree-based models such as LightGBM, Extra Trees, XGBoost, and histogram-based gradient boosting (HistGBR), demonstrated instances where 3D configurations surpassed 2D configurations in certain density classes; however, these advantages were not consistent across densities or at the overall scale. Deep learning models, including multilayer perceptron (MLP), Transformer, and graph convolutional network (GCN), exhibited a decline in performance under 3D configurations in specific density classes, and the superiority of the 3D setting was equally unstable. Overall, considering the current data scale, variable set, and modeling parameters, the combination of RF with the 3D competition-pair configuration most effectively captured the predictive improvements associated with DEM-constrained three-dimensional competition indices across density classes and at the overall scale.

3.4 Feature importance and model interpretability

After comparison of the eight candidate models, Random Forest was retained for model interpretation because it demonstrated the most stable and accurate performance across the evaluated stand-density classes. Consequently, Figures 7–10 present SHAP results for RF only. SHAP-based feature attribution was selected because it provides additive, sample-level estimates of both the direction and magnitude of predictor contributions to model output, whereas the conventional impurity-based feature importance produced by RF summarizes only the overall reduction in node impurity and does not indicate the direction of predictor contributions.

The RF-based SHAP results indicated that focal-tree growth attributes collectively made the largest contributions across the four input configurations and stand-density classes. In the 2D competition-pair configuration (Figure 7), focal-tree DBH ranked first in low-, medium-, and high-density stands as well as in the overall dataset. Additionally, focal-tree H also remained among the most important variables. Although focal-tree CW showed a relatively lower importance when considered individually, the combined mean absolute SHAP values of focal-tree DBH, H, and CW ranged from 82.2 to 92.5% across stand-density classes and reached 91.4% for the overall dataset. These dominant SHAP contributions were expected because the reference individual-tree volume was calculated from field-measured DBH and H using Equation 1, while their corresponding LiDAR-derived values were used as predictors; DBH consistently ranked first, consistent with its larger exponent in Equation 1, whereas neighbor-tree attributes and competition indices provided complementary information for individual-tree volume prediction beyond DBH and H.

Under the 2D input settings (Figures 7, 8), the Hegyi index, Arney index, and CCI showed relatively low SHAP values across most stand-density classes. However, after incorporating true 3D distance and slope weighting (Figures 9, 10), changes in the SHAP contributions of the competition indices varied among indices, stand-density classes, and input configurations. Under the competition-pair configuration (Figures 7, 9), the mean relative SHAP contributions of CCI and Arney, averaged across the three stand-density classes, increased from 2.63 to 3.67% and from 1.13 to 1.73%, respectively, whereas that of Hegyi remained at 1.23% because increases in low- and medium-density stands were offset by a decrease in high-density stands. In medium- and high-density stands, some three-dimensional competition indices approached or exceeded individual neighbor variables, while focal-tree DBH and H remained dominant. These normalized SHAP values represent relative contributions to RF prediction.

4 Discussion

The results indicate that DEM-constrained three-dimensional competition indices can provide additional neighborhood information, although the magnitude of improvement depended on input organization, stand-density class, and model choice.

4.1 Data fusion and single-tree segmentation

The integration of BLS and UAV LiDAR yielded generally stable geometric fusion across stand densities and canopy conditions (Shao et al., 2022; Zhang et al., 2025). Although stronger stem occlusion and canopy overlap in high-density stands may make point-cloud registration and stem sampling more challenging, the centimeter-level residuals observed across the plots suggest that these conditions did not preclude reliable single-tree structural reconstruction. Incomplete stem sampling caused by understory occlusion and canopy overlap may still introduce uncertainty into stem-geometry retrieval and the competition indices derived from neighbor detection and inter-tree distances. The high recall and precision of individual-tree segmentation indicate that stems were reliably detected across plots, providing a solid basis for subsequent attribute retrieval. Nevertheless, the remaining omission and commission errors may, respectively, exclude true neighboring trees or introduce spurious neighbors, thereby affecting neighbor identification, neighbor counts, inter-tree distances, and the resulting competition indices. LiDAR-derived DBH, H, and CW showed good agreement with field measurements, suggesting that the fused point clouds preserved both stem geometry and crown structure (Li Y. et al., 2025). Importantly, structural-attribute retrieval accuracy did not vary monotonically with stand density. The higher pooled R² for DBH likely reflects the dense lateral stem information provided by BLS and the clearly defined measurement height of 1.3 m. Recent comparisons of under-canopy laser-scanning systems have likewise reported high accuracy for DBH and stem-curve retrieval, including under difficult forest conditions (Muhojoki et al., 2024). In contrast, H and CW retrieval depended on the complete representation of treetops and crown boundaries. Occlusion can result in incomplete capture of upper-tree components, while contact and interlocking between adjacent crowns complicate crown delineation (Bienert et al., 2021; Tao et al., 2015). Complex stand structures characterized by dense canopy and abundant understory vegetation can also reduce tree-parameter extraction performance in Chinese fir plantations (Li Y. et al., 2025). Differences in field measurement procedures may have further contributed to the observed accuracy patterns. Field DBH measurement with a diameter tape at a fixed height was relatively direct. By comparison, field tree-height measurements can be affected by difficulty in identifying the actual treetop, whereas crown dimensions derived from simple field measurements can be biased by irregular crown geometry and simplifying assumptions about crown shape (Bienert et al., 2021; Owen and Lines, 2024). Therefore, uncertainty in the field reference measurements, together with incomplete point-cloud representation of treetops and crown boundaries, may partly explain the lower agreement observed for H and CW. A UAV-only baseline was not included in the present analysis, so the incremental contribution of BLS to attribute retrieval cannot be quantified directly and should be tested explicitly in future work. Overall, the fused point clouds supported robust parameter extraction in structurally complex stands, enabling subsequent competition-index calculation and volume modeling.

4.2 Differences between 2D and 3D competition indices

The comparison of 2D and 3D competition indices revealed systematic differences in how competitive relationships were represented under the two input organizations (Rieder et al., 2024). After incorporating DEM-derived elevation differences and slope-related weighting, the 3D indices generally showed lower median values than the corresponding 2D indices, particularly under the competition-pair configuration. Under the structural-unit configuration, the single-neighbor contribution terms within the Hegyi, Arney, and CCI summations were averaged across all neighbors of each focal tree, yielding one mean contribution value per formulation per focal tree and consequently producing more compact distributions. These lower values indicate that, after incorporating elevation differences and slope-related weighting, the same focal–neighbor relationships were represented with less influence from horizontal proximity than in the 2D indices. Across stand-density classes, high-density stands generally showed higher index values, while the density-class ordering varied among indices and input configurations. Overall, these results indicate that DEM-constrained three-dimensional competition indices modify the distributional characteristics of competition indices and provide an alternative basis for subsequent modeling of neighborhood competition (Rieder et al., 2024).

4.3 Model performance differences and the relative stability of RF

Model performance was sensitive to both input dimensionality and how neighborhood competition was organized. Overall, replacing conventional two-dimensional competition indices with DEM-constrained three-dimensional competition indices produced the clearest gains under the competition-pair configuration, particularly for RF. For RF, R² increased from 0.74 to 0.85 and RMSE decreased from 0.0434 to 0.0346 m³ tree⁻¹ at the overall scale; in high-density stands, R² increased from 0.81 to 0.90 and RMSE decreased from 0.0465 to 0.0266 m³ tree⁻¹. This suggests that in dense stands, explicit 3D distance and slope information can be effectively translated into improved single-tree volume prediction (Long et al., 2020; Rieder et al., 2024). Together with the SHAP results, this pattern suggests that terrain-aware neighborhood variables provided complementary predictive information beyond focal-tree size attributes in the examined medium- and high-density stands. However, because the changes in relative SHAP contributions varied among indices and input configurations, this pattern should not be interpreted as evidence of a uniform increase in predictor informativeness with stand crowding or as direct evidence of stronger biological competition. The more consistent gains of the competition-pair configuration may reflect its preservation of neighbor-specific variability that was partly reduced by averaging in the structural-unit configuration. The overall RF performance under the 3D competition-pair configuration also exceeded the tree-level UAV LiDAR results reported for Chinese fir by Yang Z. et al. (2025; R² = 0.639–0.725; RMSE = 0.050–0.058 m³), representing an R² increase of 0.125–0.211 and an RMSE reduction of 30.8–40.3%. Previous studies have integrated handheld and airborne LiDAR to retrieve Chinese fir tree parameters and stem volume. Other studies have represented neighborhood competition using inventory-based or 3D point-cloud-derived indices, as well as object-based or point-cloud-based competition models. However, these studies did not conduct the controlled comparison between conventional 2D and DEM-constrained 3D competition indices used in the present study (Rieder et al., 2024; Ronoud et al., 2024; Li Y. et al., 2025). By contrast, LightGBM, Extra Trees, XGBoost, and HistGBR showed less consistent advantages across density classes and input organizations, and some configurations exhibited reduced R² under the 3D structural-unit setting. Deep-learning models, including the MLP, Transformer, and GCN, were even more sensitive to input organization and showed pronounced variability under 3D inputs. Compared with RF, the other tree-based and deep-learning models showed less consistent changes across input configurations. This pattern is broadly consistent with tabular-data benchmarks showing that tree-based methods often remain competitive with neural networks, although their relative performance depends on dataset characteristics and model tuning (Uddin and Lu, 2024). Nevertheless, the coupled covariation between diverse tree growth variables and multifaceted competition indices serves as a key factor inducing cross-fold prediction fluctuations and occasional inferior performance of high-capacity deep learning models compared with RF (McElfresh et al., 2023). RF may have achieved this robustness by using dominant focal-tree size predictors (DBH, H, and CW) as a strong baseline while selectively incorporating complementary information from the 3D competition indices. The limited net change in the DBH-based Hegyi contribution likely reflects the dominant role of DBH-related information already represented in the RF model, whereas the increased relative contributions of the 3D CCI and Arney indices suggest that terrain-adjusted crown-width- and height-based neighborhood information provided additional predictive information under the competition-pair configuration. Overall, these results highlight that improving volume prediction with 3D structural information requires not only informative DEM-constrained three-dimensional competition indices but also an input organization and model choice that can reliably convert three-dimensional structural information into stable performance gains.

4.4 Drivers of Forest stock volume and the significance of terrain

Forest stock volume reflects the combined influence of environmental conditions, stand structure, and forest management practices operating across spatial scales. At broader spatial scales, Chinese fir growth responds to climatic variation, while stand density and climate jointly influence carbon storage in Chinese fir plantations (Jiang et al., 2022; Li H. et al., 2025). At local scales, topographic variation can alter the geometric representation of focal-tree–neighbor relationships and is associated with spatial differences in soil microclimate, soil moisture, soil nutrients, and understory vegetation (Nie et al., 2019; Fan et al., 2020; Liu Y. et al., 2020). Previous research has also demonstrated that incorporating slope, aspect, and competition information can improve tree height–diameter modeling, illustrating the potential value of combining terrain and neighborhood information in individual-tree prediction (Long et al., 2020).

In the present study, the role of terrain was evaluated by incorporating DEM-derived elevation differences, three-dimensional inter-tree distances, and slope weighting into the representation of focal-tree–neighbor relationships. Replacing conventional two-dimensional competition indices with their DEM-constrained three-dimensional counterparts produced the clearest and most consistent improvement in RF prediction under the competition-pair configuration. Across the slope conditions represented by the three stand groups, the magnitude of the 2D-to-3D improvement did not simply increase with mean plot slope. The low-, medium-, and high-density groups had mean slopes of 24.3°, 22.9°, and 28.1°, respectively, while the corresponding increases in R² under the RF competition-pair configuration were 0.14, 0.20, and 0.09. This pattern is consistent with the formulation of the DEM-constrained 3D indices, because terrain is represented at the focal-tree–neighbor level through elevation differences, true three-dimensional distances, and slope weighting rather than by plot-level mean slope alone. In this way, the 3D formulation refines the geometric representation of local neighborhood relationships on uneven terrain, consistent with previous findings that incorporating terrain and competition information can improve individual-tree prediction (Long et al., 2020). Together with the SHAP results, this improvement indicates that terrain-aware neighborhood variables provided complementary predictive information beyond focal-tree size attributes, particularly in the examined medium- and high-density stands. Overall, these results highlight the value of terrain-aware neighborhood representation for fine-scale individual-tree volume estimation and provide a basis for future studies combining three-dimensional structural information with direct observations of resource use and tree competition to further investigate the ecological mechanisms underlying terrain-related neighborhood relationships.

5 Conclusion

This study used fused UAV LiDAR and BLS point clouds from 12 Chinese fir plantation plots spanning contrasting density classes in Fengshushan Forest Farm, Jiangxi Province, China, to extract single-tree structural attributes and construct DEM-constrained three-dimensional competition indices for individual-tree volume prediction. (1) Compared with their corresponding two-dimensional indices, the three-dimensional indices generally showed lower median values, while the ordering across stand-density classes varied among indices and input configurations. (2) Under the competition-pair configuration, replacing the conventional two-dimensional competition indices with their DEM-constrained three-dimensional counterparts improved Random Forest performance across all stand-density classes, with the largest R² increase observed in medium-density stands and the largest RMSE reduction observed in high-density stands. Among all evaluated model–input configurations, Random Forest with the three-dimensional competition-pair input produced the clearest and most consistent predictive improvement. (3) SHAP analysis indicated that focal-tree size attributes made the largest contributions to individual-tree volume predictions, with DBH and height contributing the most and crown width providing additional predictive information. Neighbor-tree attributes generally formed a secondary group of predictors, while the three-dimensional competition indices provided clearer complementary predictive information in medium- and high-density stands. These findings indicate that pairwise focal-tree–neighbor information incorporating DEM-derived elevation differences, three-dimensional inter-tree distances, and slope weighting can complement focal-tree size attributes in fine-scale individual-tree volume estimation. Future studies could further validate the applicability of the variable set used in this study for individual-tree volume estimation in Chinese fir plantations across different cultivation regions and in plantations of other tree species. Variables derived from spaceborne LiDAR and spaceborne radar could also be incorporated to explore the potential application of this approach at a broader spatial scale.

Generative AI statement

The author(s) declared that Generative AI was not used in the creation of this manuscript.

Any alternative text (alt text) provided alongside figures in this article has been generated by Frontiers with the support of artificial intelligence and reasonable efforts have been made to ensure accuracy, including review by the authors wherever possible. If you identify any issues, please contact us.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Footnotes

Edited by: Yuhui Weng, Stephen F. Austin State University, United States

Reviewed by: Jiyeon Kim, Stephen F. Austin State University, United States

Zilong Qin, Taiyuan University of Science and Technology, China

Chongzhi Chen, Zhejiang Agriculture and Forestry University, China