Spatial Effects Analysis on Individual-Tree Aboveground Biomass in a Tropical Pinus kesiya var. langbianensis Natural Forest in Yunnan, Southwestern China

Abstract

It is essential to analyze the spatial autocorrelation and heterogeneity of aboveground biomass (AGB). But it is difficult to accurately describe due to the lack of data in clear-cutting plots. Thus, measuring the AGB directly in a clear-cutting plot can provide a reference for accurately describing the spatial variation. Therefore, a 0.3-hectare clear-cutting sample plot of Pinus kesiya var. langbianensis natural forest was selected, and the AGB was calculated by each component. The intra-group variance was quantitatively described in terms of spatial heterogeneity, and the spatial autocorrelation was explored by global and local Moran’s I. The results indicated that (1) there was different spatial heterogeneity for the different trees and organs. The intra-group variance tended to be stable after 20 m for P. kesiya var. langbianensis (PK) and other upper trees (UPs) and after 10 m for the other lower trees (LTs). (2) The spatial autocorrelation of AGB and wood biomass was similar, while the bark biomass and foliage biomass were consistent. PK and other UPs also exhibited strong spatial autocorrelation, with maximum Moran’s I values of 0.1537 and 0.1644, respectively. (3) There was spatial heterogeneity in the different components except for the bark of PK. The lowest spatial heterogeneity was found for LT.

1. Introduction

Many factors affect forest growth, such as competition with other trees, tree distribution, forest disturbance, etc. [1,2,3]. The interaction among these trees produces spatial effects. Forestry data generally have spatial effects, making them relevant and independent. Without considering spatial effects, it may lead to deviations or errors in the analysis of biomass change rules and estimation [4,5,6]. Generally, the spatial effects are usually described in terms of spatial autocorrelation and spatial heterogeneity [7,8].

Spatial autocorrelation analysis is one of the most commonly used methods for spatial analysis, such as spatial structure, clustering test analysis, recognition of spatial distribution characteristics, and selection of appropriate spatial scales [9]. Generally, it can be divided into global and local spatial autocorrelation [10]. Global spatial autocorrelation mainly describes the average level of spatial objects over the whole research area. Moran index [11], Geary’s C index [12], or Getis’G statistics [13] are commonly used to detect global spatial autocorrelation [8,14,15]. Local spatial autocorrelation can test whether there is similar or different observation aggregation in the local area, then determine the specific location of the hot spot, cold point, and outlier. The local indicators of spatial association (LISA), Getis’G statistics, and Moran scatter plots are the main indices for local spatial autocorrelation analysis [16,17]. These have been widely applied to improve the understanding of forest distribution and estimate the accuracy of forest information in forestry [18,19,20,21,22,23,24].

Spatial heterogeneity refers to the presence of structural instability characterized by systematic variations in model parameters or different response functions within a spatial dataset [4,7]. Neglecting spatial heterogeneity can result in erroneous significance tests and suboptimal predictions [16]. The importance of spatial heterogeneity lies in its provision of crucial information and insights in various fields, such as ecology, resource management, disaster risk assessment, and urban planning. It aids in better understanding and addressing spatial changes and patterns. For instance, effective forest resource management requires considering the spatial heterogeneity of biomass and growth rates to develop appropriate forestry management plans [25]. For forestry data, it represents spatial complexity, variation of ecosystem characteristics, and spatial non-homogeneity or the differences caused by the non-random distribution of forest vegetation [26,27,28,29,30,31,32,33]. Fully understanding and utilizing spatial heterogeneity will help us better understand the forest vegetation growth process and even the forest ecosystem’s evolution [34,35,36,37]. For different spatial data, the methods to evaluate and measure spatial heterogeneity include a random index, aggregated index, the nearest distance, trend surface, spectrum analysis, variation function, fractal component, analysis of variance ratio, correlation analysis, etc. [38,39,40,41,42,43,44]. Moreover, spatial autocorrelation and spatial regression models could also be used to analyze spatial heterogeneity [11,12,13,16]. Furthermore, spatial regressions, particularly the geographically weighted regression model (GWR) model, are still the most common method to reflect spatial heterogeneity [7,45,46,47,48].

Biomass is an essential attribute of forests, and there are some studies on spatial effects analysis and the spatial regression model construction of forest biomass at different scales. The application of remote sensing, geographic information systems, and global position system technologies make the spatial location information of data relatively easy to obtain and the spatial effect analysis at the regional level [49]. At the stand scale, spatial regression models, considering spatial effects, can also play a role in fitting related variables at the stand scale [50]. Due to the spatial effects of forest data, forests are influenced by their interactions with the surrounding environment during their growth process. Therefore, it is assumed that spatial regression models can meet the requirements of data independence [51]. The AGB of a forest was estimated using an ordinary least squares model (OLS), spatial lag model (SLM), spatial error model (SEM), and GWR model, and the GWR model has good performance and can solve the problems of high underestimation and low overestimation to some extent [52,53]. At the individual tree scale, Ou et al. [54] constructed the individual tree biomass estimation model of Pinus kesiya var. langbianensis by geographic weighting regression, and the fitting accuracy was significantly higher than that of the ordinary regression model. Some researchers also constructed the mixed-effect model of individual tree biomass with higher fitting and prediction accuracy by using the mixed-effect model technology. These mixed-effect models do not consider the location information of their data but mostly use some functions to represent spatial autocorrelation in the selection of variance and covariance matrix within the group [5,55]. It can be found that spatial effects are common in forest biomass data. Considering the spatial effects of individual tree biomass, it can more accurately describe the change law of biomass. Additionally, it can improve the accuracy of forest biomass estimation to a certain extent.

Regarding the spatial effect analysis of the growth of individual trees, most studies focus on the spatial effect analysis of individual tree diameters and the relationship between the diameter and tree height. For example, Zhang et al. [5] analyzed the multivariate relationship of spatial heterogeneity between tree growth and diameter; Lu and Zhang [56,57,58] used a spatial regression model to fit tree height and diameter growth under the influence of spatial effect. Gu and Zhao [59] established a GWR model for forest growth and achieved high fitting accuracy. It can be seen that the individual tree diameter at breast height (DBH) and tree height (H) has been confirmed to have spatial effects, and most of these studies focus on the changes in diameter and tree height related to volume [60,61,62]. Moreover, the forest structural heterogeneity and spatial patterns of trees associated with age structure have been reported previously [63,64], but spatial effects on individual tree biomass have rarely been reported.

It can be seen that the spatial effects of biomass and carbon storage at a regional scale and the spatial effects of DBH and the volume of trees at the individual or stand levels have been analyzed. Previous studies have focused on the change in diameter and height related to volume, which can reflect the corresponding wood or wood biomass to a certain extent. However, a tree comprises wood, bark, branches, foliage, roots, and other components. Pearce et al. [65] used a linear mixed-effect model and site as a random effect to estimate the biomass of a New Zealand shrubbery firewood forest. Chai [66] proposed a biomass ratio to describe the spatial structure of typical secondary forests in the mid-mountain zone of the Qinling Mountains based on the relationship between adjacent trees. However, the biomass data used in the quantitative evaluation were not measured but calculated using the allometric growth equation of biomass of predecessors. Unfortunately, it is difficult to obtain accurate measurement data for forest biomass, and most of the biomass data investigated in the past do not record spatial information. Therefore, the variation in spatial effects of each component regarding individual tree biomass within a stand and the differences in spatial effects of individual tree biomass among different stand types has rarely been reported. Therefore, the differences in the spatial effects among the components need to be further studied.

Then, the spatial heterogeneity and autocorrelation were analyzed at an individual tree scale. The significant contributions of this work are:

(1) To describe the spatial heterogeneity and autocorrelation of AGB at an individual tree scale using a clear-cutting plot;

(2) To explain the differences in spatial effects among the different components of the AGB;

(3) To explain the differences in spatial effects among the different tree species in a tropical natural forest.

2. Materials and Methods

2.1. Study Site

The study area is located in the town of Yutang in Mojiang County, which belongs to Pu’er City, Yunnan Province, in southwest China (Figure 1). Pinus kesiya var. langbianensis forest is one of the area’s main forest types [67]. Mojiang County is located between 101°08′–102°04′ E and 22°51′–23°59′ N. The lowest altitude is 440 meters, and the highest altitude is 2278 meters. It has a typical subtropical monsoon climate, and the seasons here are not so obvious, but there is a distinct rainy and dry season. In addition, there was abundant rainfall, with an average annual rainfall of 1338.0 mm. The rainy season happens from May to October. According to the meteorological data (1980–2010) of Mojiang station, the mean annual temperature was about 17.8 °C. The coldest month was January, with an average temperature of 11.7 °C, and the hottest month was June, with an average temperature of 21.9 °C. Soil type is mainly red soil or lateritic red soil [68].

Pinus kesiya var. langbianensis is a tree species of the Pinus genus belonging to the Pinaceae family, and it is a geographic variation of Pinus kesiya [68,69]. The Pinus kesiya var. langbianensis forests, one of the specific types in Yunnan Province, were mainly distributed in the southern subtropical zone of the west Ailao mountains in Yunnan Province. Its stand volume accounted for about 11% of the forest land of Yunnan Province. Moreover, it is also distributed in central and northern Vietnam and Laos, etc. [70]. Rapid growth is a critical feature because its trunk and annual branches grow two or more rounds per year compared to other coniferous species. Therefore, it has been an important afforestation species in Yunnan [67].

2.2. Data Investigation

A typical sample plot of 0.3 hm² (30 m × 100 m) was investigated in April and May 2015. Co-ordinate information on the location, elevation, slope degree, and slope aspect of the sample plot was recorded, and the center co-ordinates are 23°09′58.2″ N, 101°29′14.4″ E (

Table 1. The basic characteristics of the clear-cutting plot.

Location	Latitude	Longitude	Altitude (m)	Slope Degree (°)	Aspect of Slope (°)
Yutang	23°9′58.2″ N	101°29′14.4″ E	1530	22	298

). In the plot, the Pinus kesiya var. langbianensis natural forest is a natural mature mixed forest. The upper layer is dominated by Pinus kesiya var. langbianensis and some hard broadleaf species, such as Castanopsis hystrix, C. indica, and Cyclobalanopsis spp., and other soft broadleaf species, such as Schima wallichii and Alnus nepalensis. The lower layer consists of Wendlandia sp., Phyllanthus emblica, Anneslea fragrans, etc. The forest floor was mostly covered with fallen foliage and branches.

2.3. Biomass Measure and Calculate

We numbered each tree starting from the southwest corner of the sample plot and recorded the position of all trees; then, all trees were felled. Then, the DBH of trees was measured using a diameter tape measure. Meanwhile, we also measured the height of the fallen sample trees. The location, height, and DBH of each tree were recorded for all live trees with DBH > 3 cm and height ≥ 2 m (Figure 2). All trees could be divided into three groups, Pinus kesiya var. langbianensis (PK), other trees at the upper layer (UP), and the trees at the lower layer (LT). There were 132 PK, 119 UP, and 261 LT. UP includes 17 Castanopsis hystrix, 17 C. indica, 16 Cyclobalanopsis spp., 54 Schima wallichii, 12 Alnus nepalensis, and three other trees. In the LT, there are 93 Wendlandia sp., 128 Phyllanthus emblica, 32 Anneslea fragrans, and eight others species.

Biomass data for tree components were obtained using the weighing method and the bulk density method, including wood, bark, branches, foliage. For wood and bark, the discs of individual trees with different heights and bark thicknesses were intercepted and weighed for fresh weight to obtain moisture content and wood density. For foliage and branches, the branches and foliage were weighed in the field and in the sample to assess dry weight [67]. All samples were dried in the laboratory at 105 °C to a fixed mass. The density of the disc and the corresponding volume of each stem segment were used to calculate the dry weight of the wood and bark using the proportion of the total fresh mass and the moisture content of the sample to determine the dry weight of the branches and foliage [71,72]. The information is listed in

Table 2. The basic information of the sampling trees in the clear-cutting plot.

Trees or Trees Group		Age (Years)	Height (m)	DBH (cm)	Biomass (kg)
					Wood	Bark	Branches	Foliage	Aboveground
Pinus kesiya var. langbianensis (n = 132)	Mean	37	17.05	24.26	236.44	6.78	37.24	4.45	284.91
	Std. Err.	10	4.58	9.93	224.73	5.33	47.83	4.74	268.62
	Min.	16	6.80	7.00	3.29	0.14	0.51	0.06	4.07
	Max.	66	25.60	45.10	1075.79	34.49	307.23	23.39	1180.42
Other upper trees (n = 119)	Mean	25	9.97	14.76	55.24	15.16	13.62	3.45	87.47
	Std. Err.	10	4.21	7.76	72.61	19.09	17.78	4.97	102.08
	Min.	5	3.30	3.40	0.82	0.29	0.07	0.02	1.74
	Max.	48	21.50	36.00	368.49	88.85	103.17	31.70	482.72
Other lower trees (n = 261)	Mean	16	6.36	6.73	5.97	1.32	1.77	0.54	9.60
	Std. Err.	5	1.66	2.62	6.81	1.55	2.74	0.89	11.09
	Min.	5	2.20	4.00	0.83	0.07	0.03	0.00	1.22
	Max.	35	12.00	21.00	44.20	10.57	21.04	7.67	71.80
Total (n = 512)	Mean	24	9.95	13.12	76.84	5.95	13.67	2.22	98.68
	Std. Err.	12	5.53	9.81	153.31	11.14	29.71	3.88	184.67
	Min.	5	2.20	3.40	0.82	0.07	0.03	0.00	1.22
	Max.	66	25.60	45.10	1075.79	88.85	307.23	31.70	1180.42

.

2.4. Spatial Effects Analysis

2.4.1. Spatial Heterogeneity of the Aboveground Biomass

The intra-group variance can be used to describe the spatial heterogeneity of study subjects quantitatively. This is a function of block size and usually increases as the size of the blocks increases [73], and it can be calculated by Equation (1).

$$S_{\mathrm{intra}}=\frac{1}{B}{\displaystyle \sum _{g=1}^B\frac{1}{n_g}}{\displaystyle \sum _{h=1}^{n_g}{\left(X_{gh}-{\overline{X}}_g\right)}^2}$$

(1)

where B is the number of groups at a particular grouping distance; X_gh is the h-th observation biomass value in the g-th block; ${\stackrel{\text{−}}{X}}_g$ is the mean of observation value in the g-th block.

2.4.2. Spatial Autocorrelation of the Aboveground Biomass

Moran’s I effectively reflect the differences and correlations in the spatial distribution of observations. It can also reflect the overall spatial clustering pattern of objects in the study area [11]. The value range of Moran’s I is between (−1, 1). The closer the absolute value of the index value is to 1, the greater the degree of aggregation. When Moran’s I is positive, the observations within the distance (h) tend to be similar.

On the contrary, they tend to be dissimilar when the coefficient is negative. If Moran’s I index is close to or equal to 0, the observations are spatially random and independent. So, Moran’s I was chosen to describe the spatial autocorrelation for each component biomass (Equation (2)).

$$I=\left(\frac{n}{{\displaystyle \sum _{i=1}^n{\displaystyle \sum _{j=1}^n{w}_{ij}\left(d\right)}}}\right)\frac{{\displaystyle \sum _{i=1}^n{\displaystyle \sum _{j=1}^n{w}_{ij}\left(d\right)}(x_i-\overline{x})(x_j-\overline{x})}}{{\displaystyle \sum _{i=1}^n{(x_i-\overline{x})}^2}}$$

(2)

where n is the number of plots; x_i and x_j are the observation value at plot i and plot j, $\overline{x}$ is the average of all attribute values; $w_{ij}\left(d\right)$ is the spatial weight matrix value that represents the relationship between plot i and plot j of the spatial object.

After the Moran index value is calculated, the significance of the index is tested by a Z-score. The higher (or lower) the Z-score, the higher the degree of clustering. A positive Z-score indicates a cluster of high values, and a negative Z-score indicates clusters of low values. If the Z-score is close to zero, no significant clustering exists within this area, and the Z-score was calculated by Equation (3).

$$Z(I)=\frac{I-E(I)}{\sqrt{Var(I)}},E(I)=\frac{-1}{n-1}$$

(3)

where Z(I) is an index to measure the strength of spatial clustering patterns; E(I) and Var(I) are the expected value and variance of the index value I, respectively;

Moreover, local spatial autocorrelation can be used to describe the degree of spatial correlation between a research object and its adjacent units. In this study, the distance corresponding to the first peak after each component’s spatial clustering pattern for biomass reached significance was used as the bandwidth, and ArcGIS was used to calculate local Moran’s I_i. Furthermore, the local Moran’s I_i was calculated by Equation (4):

$$I_i=\frac{n^2}{{\displaystyle \sum _i^n{\displaystyle \sum _j^n{w}_{ij}}}}\frac{(x_i-\overline{x}){\displaystyle \sum _j^n{w}_{ij}(x_j-\overline{x})}}{{\displaystyle \sum _i^n{(x_i-\overline{x})}^2}}$$

(4)

where n is the number of plots; x_i and x_j are the observation value at plot i and plot j, $\stackrel{\text{−}}{x}$ is the average of all attribute values; $w_{ij}\left(\mathrm{d}\right)$ is the spatial weight matrix value. Different from the global Moran’s I, the value range of the local Moran’s I_i is not limited to (−1, 1); if the I_i value is positive, it shows that the location of the plot is a positive correlation and reflects the aggregation of similar values. But the I_i is negative, and it shows that the location of the plot is a negative correlation and reflects the aggregation of different values.

The global Moran’s I and Z (I) values were calculated using the incremental spatial autocorrelation tool in ArcGIS. The most significant distance of the spatial clustering pattern obtained from the incremental spatial autocorrelation analysis was used as the scale parameter for local spatial autocorrelation analysis. On this basis, using ArcGIS clustering and outlier analysis tools, the spatial weight matrix was constructed in the form of fixed distance bandwidth, and the local Moran’s I_i was calculated to quantitatively describe the spatial correlation between local, regional observations.

3. Results

3.1. Spatial Heterogeneity of the Aboveground Biomass

The spatial heterogeneity variation in biomass of each component showed the same trend in the range of 0–100 m (Figure 3). Aboveground, wood, branches, bark, and foliage were the order of the intra-group variance values for the different components. Although the value changed with the different organs, the spatial variation gradually increased with increasing distance; after reaching a certain distance, the variation tended to be stable. The variance values tended to be stable after 20 m for both PK and UP. In comparison, the variation of LT and all tree species tended to be stable after 10 m. Furthermore, the spatial variance values of PK were significantly larger than the other tree types for the AGB, wood, and branches. The spatial variance values of UP were the largest for both bark and foliage, and all the biomass components of the LT showed the smallest spatial variance values.

3.2. Spatial Autocorrelation of the Aboveground Biomass

3.2.1. Global Moran’s I

The global Moran’s I of each component biomass at the different lag distances were calculated, and the Moran’s I correlogram and Z-score values of the biomass with lag distances are shown, respectively. For all trees (Figure 4), the distances at the first peak (DFP) with a significant Z-score larger than 1.96 were 21.2 m for the AGB and wood biomass. The distance value was 6.2 m, and the biomass for both the foliage and branch was 4.8 m. Moreover, at the corresponding lag distance, the spatial autocorrelation was significant because the Z-scores were larger than 1.96. Meanwhile, the Moran’s I values of AGB and wood were lower than 0.01, and the values of foliage, bark, and branches were larger than 0.034.

For the different tree types, the PK had a similar change law for all trees; the DFP was 12.4 m for both AGB and wood biomass, 10.0 m for bark biomass, and 7.4 m for both foliage and branches. Furthermore, the Moran’s I was larger than 0.13 for foliage and branches; the values for the AGB, wood, and bark were 0.0687, 0.0609, and 0.0389, respectively (Figure 5). For the UP, the spatial autocorrelation of both foliage and branches was not significant because the Z score was not significant along with the lag distance, and the absolute Z scores were not larger than 1.96. Meanwhile, the DFP of AGB, wood, and bark were 6.6 m, 6.6 m, and 9.0 m, respectively. The corresponding Moran’s I values were 0.1277, 0.1372, and 0.1644 (Figure 6). Moreover, for LT, the branch biomass had no significant spatial autocorrelation, and the Moran’s I of both AGB and wood were 0.0520 and 0.0617 at a DFP of 6.4 m. The bark biomass had the largest DFP at 17.5 m, and the Moran’s I was 0.0233. Meanwhile, the foliage biomass had the lowest DFP at 6.0 m, and the Moran’s I was 0.0771 (Figure 7).

Moreover, we can find that the AGB and wood biomass had the same DFP for each tree type for all trees. Additionally, they had the highest DFP for all trees and PK but the lowest UP, and most of the foliage biomass had a lower DFP, with the largest Moran’s I for each tree type and all trees. Furthermore, a large portion of the PK and UP had the highest Moran’s I when compared to that of LT.

3.2.2. Local Moran’s I_i

The counting statistics of spatial autocorrelation were made and listed in

Table 3. Spatial distribution by Z-score using all trees. PK is Pinus kesiya var. langbianensis, UP is the other trees at the upper layer, and LT represents the other trees at the lower layer. NS is the point with no significant spatial autocorrelation, HH is the similar aggregation point with high-high clustering, LL is the similar aggregation point with low-low clustering, HL is the dissimilarity aggregation point with a high-low outlier, and LH is the dissimilarity aggregation point with a low-high outlier. N is the number of degrees of distribution for each component, and P is the proportion of distribution of each component.

Tree Types	Components	NS		HH		HL		LH		LL
		N	P (%)	N	P (%)	N	P (%)	N	P (%)	N	P (%)
PK	Wood	102	77.27	11	8.33	19	14.39	0	0.00	0	0.00
	Bark	128	96.97	4	3.03	0	0.00	0	0.00	0	0.00
	Branches	118	89.39	8	6.06	5	3.79	1	0.76	0	0.00
	Foliage	116	87.88	15	11.36	1	0.76	0	0.00	0	0.00
	Aboveground	103	78.03	12	9.09	17	12.88	0	0.00	0	0.00
UP	Wood	117	98.32	2	1.68	0	0.00	0	0.00	0	0.00
	Bark	99	83.19	15	12.61	5	4.20	0	0.00	0	0.00
	Branches	119	100.00	0	0.00	0	0.00	0	0.00	0	0.00
	Foliage	115	96.64	1	0.84	3	2.52	0	0.00	0	0.00
	Aboveground	117	98.32	2	1.68	0	0.00	0	0.00	0	0.00
LT	Wood	261	100.00	0	0.00	0	0.00	0	0.00	0	0.00
	Bark	258	98.85	0	0.00	0	0.00	3	1.15	0	0.00
	Branches	259	99.23	0	0.00	0	0.00	2	0.77	0	0.00
	Foliage	260	99.62	0	0.00	0	0.00	1	0.38	0	0.00
	Aboveground	261	100.00	0	0.00	0	0.00	0	0.00	0	0.00
ALL	Wood	480	93.75	13	2.54	19	3.71	0	0.00	0	0.00
	Bark	485	94.73	19	3.71	5	0.98	3	0.59	0	0.00
	Branches	496	96.88	8	1.56	5	0.98	3	0.59	0	0.00
	Foliage	491	95.90	16	3.13	4	0.78	1	0.20	0	0.00
	Aboveground	481	93.95	14	2.73	17	3.32	0	0.00	0	0.00

. The PK had a lower percentage of the point with no significant spatial autocorrelation (NS) when compared with UP and LT. Especially for the LT, the percentage of the NS point was greater than 98%, and all points were the NS for the wood biomass and AGB. But for the PK, the percentage with NS for the wood biomass and AGB were the lowest, and the percentage was 77.27% and 78.03%, respectively. Moreover, for the point with significant spatial autocorrelation, there was no similar aggregation point for low-low clustering (LL) for all trees; there were similar aggregation points for high-high clustering (HH), and a dissimilarity aggregation point for high-low outliers (HLs) were not found in the LT. The point with significant spatial autocorrelation in the plot was primarily LH outliers and HH clustering, and wood and AGB showed similar local spatial autocorrelation. HH points were mainly found in the wood, foliage, and AGB of PK, and bark biomass of UP. Meanwhile, HL points were mainly found in the wood and AGB of PK. For the dissimilarity aggregation points with a low-high outlier (LH), the percentage was lower; all UP was not found at the point with LH, and only 1 point for branches biomass for PK was found for LH.

Moreover, the spatial autocorrelation for each biomass component of the different tree or tree groups was shown in the local Moran’s I_i scatters. For the PK, the biomass of AGB, wood, and branches mainly showed significant HH clustering on the right side of the plot, while the spatial clustering pattern of biomass in the middle and left sides showed HL outliers. In addition, the biomass of bark and foliage mainly showed HH clustering within the plots. The biomass of some components was also accompanied by the LH outlier (Figure 8). For the UP, AGB and wood biomass showed HH clustering on the right side of the plot and NS on the left side. Similarly, foliage biomass also showed HH clustering on the right side of the plot, but an HL outlier was also present on the right and left sides. In contrast, bark biomass showed mainly HH clustering in the plot; HL outliers were found on the right center and near the left edge of the plot (Figure 9). But for the LT, the components mainly showed no significant aggregation (NS) and sporadic LH outliers in the plot (Figure 10).

Furthermore, for different tree species and components, we calculated the local Moran’s I_i of each tree using all trees. The statistics are listed in

Table 4. Local Moran’s I_i for the biomass of the different components.

Tree Species	Indices	Wood	Bark	Foliage	Branches	Aboveground
Pinus kesiya var. langbianensis	Mean	−0.45	0.10	0.39	0.17	−0.42
	Std. Err.	0.23	0.07	0.20	0.21	0.23
	Skewness	0.01	2.52	2.87	3.25	−0.03
	Kurtosis	2.72	13.59	11.83	15.04	2.18
	Min.	−10.32	−2.06	−3.87	−4.58	−9.67
	Max.	8.43	4.94	12.73	14.35	8.12
Other upper trees	Mean	0.15	0.82	−0.07	−0.04	0.07
	Std. Err.	0.09	0.31	0.11	0.06	0.10
	Skewness	1.78	3.81	0.16	−2.42	0.50
	Kurtosis	12.72	17.68	9.87	10.13	4.56
	Min.	−3.02	−3.71	−6.06	−3.76	−3.06
	Max.	6.49	20.07	6.14	1.26	5.40
Other lower trees	Mean	0.26	0.08	0.09	0.11	0.28
	Std. Err.	0.05	0.05	0.05	0.04	0.05
	Skewness	−0.72	−2.05	−2.59	−4.19	−0.86
	Kurtosis	0.18	7.46	11.12	28.11	0.79
	Min.	−2.60	−4.64	−4.23	−5.00	−3.24
	Max.	1.53	1.25	1.44	1.00	1.70
Total	Mean	0.05	0.26	0.13	0.09	0.05
	Std. Err.	0.07	0.08	0.06	0.06	0.07
	Skewness	−0.56	6.64	3.15	4.67	−0.63
	Kurtosis	9.76	64.72	26.56	46.34	7.90
	Min.	−10.32	−4.64	−6.06	−5.00	−9.67
	Max.	8.43	20.07	12.73	14.35	8.12

. According to the local autocorrelation analysis of biomass in different tree species and different components, it can be seen that, for the mean local autocorrelation value of wood biomass, the autocorrelation coefficient of PK showed a negative value, and the value (−0.45) was smaller than other trees (upper trees, lower trees, and total), with a negative autocorrelation. The minimum value (−10.32) and maximum value (8.43) were consistent with the total, and the minimum value was smaller than the other upper and lower trees, while the maximum value was larger than them. However, the standard error (0.23) was the largest among all tree species. For the local autocorrelation analysis of bark biomass, it can be seen that all tree species showed a positive correlation, and the highest mean correlation value (0.82) was found in other upper trees, and the highest value (20.07) was the same as that of the total. For the mean autocorrelation of foliage and branches biomass, except for the negative correlation of other upper trees, all showed a positive correlation, and the correlation of PK was the strongest; for foliage (0.39) and branches (0.17), the maximum value was the same as the total: foliage (12.73) and branches (14.35) and was higher than the other upper and lower trees. The mean, maximum, and minimum values of the local autocorrelation coefficients of AGB exhibited the same trend as wood biomass.

4. Discussion

4.1. Spatial Heterogeneity and Autocorrelation in the Tropical Natural Forest

The forest community has high spatial heterogeneity [74]. Under heterogeneous conditions, community differentiation led to self-thinning, weakened inter-specific competition, and enhanced individual stress resistance, which was beneficial for the stability of various groups in the community and the overall development of the community. In this study, in the range of 0–100 m, the change in spatial heterogeneity tended to be stable after a distance of 20 m. The spatial distribution of wood, bark, branches, foliage, and AGB in the natural forests of Pinus kesiya var. langbianensiss was not random but was positively correlated and showed some spatial aggregation. The semi-variance model structure could indicate the spatial heterogeneity of AGB [75]. As the distance increases, the correlation between the spatial variables and their surrounding values decreases [62], and this pattern was similar to the biomass in the broadleaved forests [76].

Moreover, the high biomass values were mainly concentrated on the right side of the sample plot, and the low values were mainly concentrated on the left and middle of the sample plot. The spatial clustering characteristics objectively existed in all biomass components [77]. Local spatial correlations were more robust for aboveground, wood, and bark biomass and mostly showed high-high aggregation and high-low outliers, with the lowest correlations for branch biomass. In addition, the findings indicated that the local spatial autocorrelation of the biomass components of lower tree species was more significant compared to that of other tree species. It could cause spatial patches within the stand by the specific dynamic processes that vary by spatial correlations among the tree species and local disturbances [64] and lead to spatial heterogeneity in the stand structure as well as within and among the tree species [78,79]. Therefore, it may be responsible for the spatial effects on individual tree biomass. Based on this, it can guide the study of the spatial structure, spatial heterogeneity, spatial autocorrelation, and biomass of Pinus kesiya var. langbianensiss natural forest in large regions. Rödig et al. [80] demonstrated how this linkage allowed for quantifying the spatial variation in forest structure caused by tree-level to regional-scale disturbances. In addition, the study of heterogeneity can also combine with remote sensing methods. Han et al. [81] estimated and explored the aboveground carbon storage (AGC) of bamboo forests by using Landsat images and five dates measured in the field. The results showed that the AGC of the moso bamboo increased gradually in the five periods, and its spatial heterogeneity had moderate autocorrelation.

The results of this study show that the AGB and the biomass of each component have spatial heterogeneity and spatial autocorrelation, which can be estimated by considering spatial regression models, such as GWR, in future studies. The GWR model has performed well in estimating biomass in many studies [51,52,53,61]. However, when estimating biomass or other tree characteristics, the applicability of the model needs to be fully considered before it can be calculated. The use of the SEM model enables a direct explanation of spatial autocorrelation. Meanwhile, the incorporation of spatial dependencies and variability in the modeling process through LMM (linear mixed model) and GWR techniques allows for improved data fitting and a more precise prediction of response variables [61].

4.2. Spatial Effect of the Different Components

The following results were obtained from the study of spatial heterogeneity regarding the biomass of different components. The components, which were ranked in terms of intra-group variance values, were aboveground, wood, branches, bark, and foliage, respectively. There was an allometric relationship for plant growth, and the plants preferentially allocated biomass to organ harvesting according to the biomass allocation theory [82]. Moreover, biomass allocation may vary with stand age, species competition, etc. [83], and therefore affect the spatial effect in the different components of biomass. Competition may have an especially meaningful impact on the biomass allocation of trees. Plant competition and herbivory interactions were crucial for reducing plant production performance and changing biomass allocation patterns [84]. The biomass allocation would prefer tree crowns with little light and tree roots with limited nutrients or water [82]. Under different light environments, there was no difference in the biomass allocation of stems and roots and leaf mass per area between invasive and native species [85]. In this study, the AGB had a much more significant spatial variation in the whole stand. This was consistent with the results of Du et al. [75] research on spatial heterogeneity and the carbon contribution of AGB of moso bamboo in Anji County, Zhejiang Province; in addition, their research also indicated that spatial autocorrelation had a determinative influence on the spatial heterogeneity and spatial distribution pattern of moso bamboo’s AGC. Therefore, the study of biomass cannot ignore the impact of spatial heterogeneity and spatial autocorrelation.

In addition, it is important to carry out a spatial autocorrelation analysis of the biomass of each component. In this study, global and local spatial autocorrelation analyses were carried out for the different components of different tree species so that the spatial autocorrelation between trees in the whole stand could be more easily seen. The global spatial autocorrelation trends for each component of Pinus kesiya var. langbianensis were consistent with those of all trees in the sample plot, while the trends for other upper trees and other lower trees differed from those of all trees in the sample plot. Zhou et al. [83] conducted biomass modeling of the moso bamboo foliage as a way to accurately estimate the dynamics of foliage biomass, thus contributing to the development of accurate simulations. The spatial autocorrelation was used by Wu et al. [86] to conduct a spatio-temporal analysis of forest stock in Sichuan Province, and a mapping exercise was conducted for global spatial autocorrelation and local spatial autocorrelation.

The local spatial autocorrelation of biomass is also essential, as the distribution of habitat factors, such as light, water, and soil, are heterogeneous in natural environments [87]. Stand volume and biomass are closely related, as habitat variation can lead to the territorial and specific structural characteristics of forest volume [88]. The local spatial autocorrelation analysis revealed distinct patterns of aggregation in the sample plots for the different components studied. Specifically, the HH, HL, and LH aggregations exhibited varying location and aggregation patterns. Stand volume was a stochastic function that related to random variables and location [86]; therefore, forest biomass estimation needs to take full account of the spatial location and distribution of trees. Spatial autocorrelation to interpret and aid model outputs can identify areas that are systematically sensitive to specific model inputs and can also guide the collection of additional data to refine model predictions [89].

4.3. Spatial Effect of the Different Tree Species

The spatial variance values of PK were significantly larger than the other tree types for the aboveground, wood, and branches biomass. The spatial variance values of UP were the largest for both bark and foliage, and the biomass components of the LT showed the smallest spatial variance values. Then, significant or highly significant spatial autocorrelations were found for the biomass of all the components except for the branches of LT and the foliage and branches for UP. The results showed that there was similar aggregation (high-high clustering or low-low clustering) or dissimilarity aggregation (high-low outlier or low-high outlier) regarding biomass in the natural stand of Pinus kesiya var. langbianensis but the correlation was not strong. Therefore, the spatial autocorrelation variation patterns of AGB and wood biomass were similar, and the patterns of foliage and branches biomass were comparatively similar (Figure 4, Figure 5, Figure 6 and Figure 7).

The species characteristics and specific dynamic process make the biomass allocation among the tree species different [64]; this then leads to the forest structure and spatial heterogeneity of different components and various tree species [78,79], which may be the reasons that affect the spatial effect of individual tree biomass. In this study, the stand was dominated by the masculine species of Pinus kesiya var. Langbianensis. Pinus kesiya var. Langbianensis is a typical fast-growing heliophilous species; the spatial autocorrelation of Pinus kesiya var. langbianensis natural forest was dominated by positive spatial autocorrelation for each component’s biomass, and there was significant positive spatial autocorrelation. Therefore, the spatial correlation was more prominent in Pinus kesiya var. langbianensis in terms of different tree species.

Moreover, the local spatial autocorrelation of each biomass component for the lower tree species was stronger than that of other tree species. The foliage biomass of the lower trees and all trees, the bark biomass of other upper trees, and the branches biomass of Pinus kesiya var. langbianensis were strongly correlated. However, it was also found that the biomass of some components had insignificant spatial autocorrelation, such as the foliage biomass of the upper tree species and the foliage biomass of the lower tree species. Subsequently, competition may be one of the main reasons.

The competition strongly affected tree productivity and reproductive capacity [84,85]. Light is one of the most critical factors controlling tree survival and growth, and the canopy formed by the upper layer of trees can limit the use of light resources by the species in the lower layer [90]. When compared to the dominant and upper-layer species, the lower-layer species had the most negligible spatial variation in the biomass of each component. The upper tree species allocated more of the obtained resources to the physiological activities of the foliage [47,61], coupled with the photosynthetic capacity of the bark being more sensitive to temperature [91]. The foliage and bark in the upper tree species got more light compared to the other types, and competition was more pronounced.

4.4. The Limits of the Study

The clear-cutting sample plot of Pinus kesiya var. langbianensis natural forest is 0.3 hectares; this study carried out spatial effect analysis on the data by calculating the variance of each component biomass and Moran’s index. It concluded that there was spatial heterogeneity in the biomass of an individual tree within a certain distance. This study explored spatial autocorrelation and spatial heterogeneity among individual trees. The study site is located in the main distribution area of Pinus kesiya var. langbianensis in Yunnan Province, namely the semi-humid areas of the northern tropical and southern subtropical regions of Yunnan Province. The stand structure is basically the same, with relatively small differences in site conditions and stand age [72]. Therefore, the sample plots were representative of the specific characteristics of the region. However, the number of sample plots could be increased in the next study; thus, spatial autocorrelation and spatial heterogeneity at the stand scale could be more easily investigated. Du et al. [75] investigated 55 sample plots with a sample size of 30 m × 30 m to study the spatial heterogeneity and spatial patterns of the above-ground biomass of moso bamboo. Meanwhile, the differences in physiological, biochemical characteristics, and stand structure between heliophilous species and shade-tolerant trees and the spatial effects of forest biomass dominated by other shade-tolerant species must be further studied.

Moreover, little research about the quantitative analysis of the spatial effects of individual tree AGB at a plot scale has been conducted because the measurement of biomass data is destructive and labor-intensive [80]. In this study, we discussed the reasons for the spatial effect of the natural forest of Pinus kesiya var. langbianensis. Moreover, only one sample plot was harvested using a clear-cutting plot; to some extent, it also reflects the organ and inter-specific differences of AGB at an individual tree scale for the Pinus kesiya var. langbianensis forest. However, there is a lack of validation or evidence from the relevant research. Therefore, in future research, the techniques employed in this study can be extended to investigate the spatial effects of biomass across a larger number of sample plots and more forest types to further confirm the spatial effect of the tree biomass in a forest at an individual tree scale. Additionally, the utilization of spatial regression models can be explored as a means to explore the spatial effect analysis for biomass.

5. Conclusions

In this study, all trees in a plot with 0.3 hectares (30 m × 100 m) for Pinus kesiya var. langbianensis natural forest were entirely cut down in Mojiang County, Yunnan Province, southwestern China. Then, the AGB of each tree was measured according to different organs, such as wood, bark, foliage, and branches. The intra-group variance was applied to describe the spatial heterogeneity of biomass in each component, and the global and local Moran indices were used to analyze the spatial autocorrelation for the three types of tree species, such as Pinus kesiya var. langbianensis (PK), the other upper trees (UPs), and the other lower trees (LTs). Therefore, these are the following conclusions.

(1) There was spatial heterogeneity for the biomass of each component in the P. kesiya var. langbianensis natural forest. The spatial heterogeneity of each biomass component showed an increase with distance and stabilization at distances of about 20 m for the UP and PK, but for the LT, the distance values were 10 m.

(2) The spatial autocorrelation differences among the organs were obvious, and some components showed a similar change trend. The spatial autocorrelation trend of AGB and wood biomass were similar, while bark biomass and foliage biomass had a similar spatial variation trend.

(3) The spatial effect was various for the different trees or tree groups. The weakest spatial heterogeneity was found in LT, while PK showed strong spatial heterogeneity in all components except bark biomass. PK and UP also showed significant spatial autocorrelation, with a global Moran’s I of 0.0389–0.1537 and 0.1277–0.1644, respectively.