Developing a Method to Estimate Above-Ground Carbon Stock of Forest Tree Species Pinus densata Using Remote Sensing and Climatic Data

Abstract

Forest above-ground carbon stock (AGCS) is one of the primary ecological evaluation indicators, so it is crucial to estimate the AGCS accurately. In this research, we added the climatic and topographic factors to the estimation process by a remote sensing approach to explore their impact and to achieve more precise estimations. We hope to develop a more accurate estimation method for AGCS based on remote sensing data and climate data. The random forest (RF) method has good robustness and wide applicability. Therefore, we modeled and predicted the AGCS by RF based on sixty field sample plots of Pinus densata pure forests in southwest China and the factors extracted from Landsat 8 OLI images (source I), Sentinel-2A images (source II), and combined Landsat 8 OLI and Sentinel-2A images (source III). We added the topographic and climatic factors to establish the AGCS estimation model and compared the results. The topographic factors contain elevation, slope, and aspect. Climatic factors contain mean annual temperature, annual precipitation, annual potential evapotranspiration, and monthly mean potential evapotranspiration. It was found that the R² and RMSE of the model based on source III were better than the R² and RMSE of the models based on source I and source II. Compared to the models based on source I and source II, the model based on source III improved R² by up to 0.08, reduced RMSE by up to 2.88 t/ha, and improved P by up to 4.29%. Among the models without adding factors, the model based on source III worked the best, with an R² of 0.87, an RMSE of 10.81 t/ha, an rRMSE of 23.19%, and a P of 79.71%. Among the models that added topographic factors, the model based on source III worked best after adding elevation, with an R² of 0.89, an RMSE of 10.01 t/ha, an rRMSE of 21.47%, and a P of 82.17%. Among the models that added climatic factors, the model that added the annual precipitation factor had the best modeling result, with an R² of 0.90, an RMSE of 9.53 t/ha, an rRMSE of 20.59%, and a P of 83.00%. The prediction result exhibited that the AGCS of the Pinus densata forest in 2021 was 9,737,487.52 t. The combination of Landsat 8 OLI and Sentinel-2A could improve the prediction accuracy of the AGCS. The addition of annual precipitation can effectively improve the accuracy of AGCS estimation. Higher resolution of climate data is needed to enhance the modeling in future work.

1. Introduction

The FAO (Food and Agriculture Organization of the United Nations) reported that the global carbon stock of biomass in live forests reached 294.535 billion t, with an average of 72.6 t/ha. This represents approximately 44.5% of the total global carbon stock in global forest ecosystems, including those in dead wood, dead leaves, and soil [1]. An accurate forest carbon stock estimation allows for a more intuitive reflection of the forest carbon potential [2]. Moreover, remote sensing can be very useful for estimating forest parameters. Therefore, it is necessary to explore and develop new methods for more accurate AGCS estimation based on remote sensing. It is conducive to facilitating the sustainable balance and stable development of the global ecological environment. The forest cover of Shangri-La [3] plays an important role in the ecological balance of northwestern Yunnan, as well as in the global ecological balance.

The Landsat series data are widely employed as a standard optical remote sensing dataset for estimating the AGCS and AGB (above-ground biomass) in forest ecosystems. In the estimation studies of the AGCS and AGB, the low resolution of Landsat data is their main drawback, while Sentinel 2 data have higher resolution as well as more band information [4], and these advantages can compensate for the shortcomings of Landsat data to some extent. However, there are still a number of studies that have only used Landsat data. For example, the use of Landsat 8 OLI remote sensing data based on mixed effects can improve the accuracy of forest AGB estimation [5]. The results of the Landsat 8 OLI-based forest AGB estimation study by Li et al. also showed that this image can be well used for forest AGB estimation [6]. Zheng et al. indicated that forest AGB can also be better estimated based on Landsat images [7]. After processing Landsat long-time series image data by a filtering algorithm, Teng et al. obtained higher AGB prediction accuracy [8]. Meanwhile, Landsat images can also be applied to estimate grassland AGB [9]. In addition to the studies mentioned above, there are many other applications of Landsat images for AGB estimation [10,11,12,13,14]. Their research objects or research methods are also different. Thus, Landsat images play an important role in the estimation of AGB.

The advantage of the Sentinel-2 images over the Landsat image is their higher spatial resolution of 10 and 20 m and temporal resolution (up to 5 days) [4]. It is feasible to apply the Sentinel-2 images to estimate the forest AGB and AGCS. The higher resolution of Sentinel-2 images tends to bring higher accuracy to the estimation of the AGB and AGCS [15]. Sentinel-2 images can provide important information in forest biomass research, such as some of the vegetation indices from Sentinel-2 [16], which can effectively enable remote sensing estimation of AGB [17]. The potential of the Sentinel-2 images in vegetation biomass studies allows them to be used in grassland biomass estimation [18,19,20] and in combination with other remote sensing images [21,22]. Combining multiple remote sensing sources may solve some difficulties of using a single source data for AGB modeling [23].

Topographic factors are commonly addressed in related studies. These studies include assessing the site quality of forest land before estimating the AGB and AGCS [24,25]. Zhao et al. also estimated the AGB based on the different site classes [26]. In our previous study, the terrain ecological niche index (TNI) was added to improve the prediction accuracy of the AGCS [27]. Additionally, AGB estimation models were built by categorizing elevation and slope [3] and utilizing topographic factors as mixed effectors [28]. Adding topographic factors allows for dynamic biomass estimation over long time series data [29]. Elevation, slope, and aspect factors were added in the estimation of the AGCS and its uncertainty analysis [30], as well as in the AGB estimation of natural forests [31]. The study by Xu et al. [32] also used topographic factors to establish the AGB estimation model. It is important to note that the specific topographic factors (such as elevation, slope, and aspect) employed in these studies vary with their specific usage and research objectives.

The AGCS and AGB are also affected by climate factors that affect forest growth [33,34], as land surface temperature can have an impact on forest carbon stock and biomass [35]. Land surface temperature is negatively correlated with forest cover, especially in the dry season [36]. Therefore, land surface temperature could naturally be applied to estimate the AGCS [37,38,39]. However, climatic factors such as temperature and precipitation have yet to be practically applied to the estimation of biomass or carbon stock in the vast majority of studies. Studies on the AGCS and climatic factors have preferred to analyze the spatial and temporal correlation between carbon stock and climatic factors [40,41,42].

Most studies on forest AGB and AGCS in Shangri-La, Yunnan Province, have primarily focused on the influence of topographic factors on AGB modeling and estimation [30,32,43,44]. A few studies have addressed climatic factors [37,45,46]. Among them, Wang et al. and Cheng et al. [45,46] used the same climatic factors, solar radiation in the growing season, growing accumulated temperature, and precipitation in the growing season, and Yin et al. [37] used the land surface temperature. Although these studies used additional factors, they only utilized the Landsat images with relatively coarse spatial resolution. In addition to this, most studies of Pinus densata in the region also have used only Landsat data [8,27,28,29,47]. In contrast, far fewer AGB or AGCS studies [44,48,49] have been conducted in the region using Sentinel-2 data.

Although some studies that have used Landsat images to estimate the AGCS have added topographic or climatic factors [24,25,27,28,29,30,31,32], they were limited by the resolution of the images used and the inadequacy of single-source remote sensing. Climatic factors were used much less often than topographic factors in studies about AGB and AGCS estimation. On the other hand, the studies using Sentinel-2 images have solved the lack of resolution of Landsat images while ignoring the effect of topographic and climatic factors [15,16,17,18,20,21,22,23].

Pinus densata is a prominent tree species in the region [3], so it should get increased attention and meticulous evaluation to facilitate more informed management decisions. Since most studies on the AGB or AGCS of Pinus densata do not consider the comprehensive influence of environmental factors and remote sensing data sources, it is evident that further supplementation is needed in the relevant research on the AGCS of the Pinus densata forest in Shangri-La. Based on this, we studied Pinus densata in Shangri-La, Yunnan Province, China. We combined the Landsat 8 OLI and Sentinel-2A images in 2021, DEM data, sample plots, mean annual temperature (MAT), annual precipitation (AP), annual potential evapotranspiration (APET), and monthly mean potential evapotranspiration (MMPET) data in 2021 of the area. However, in the Landsat images, there is a difference in the acquisition date between “LC81310412021082LGN00” and the rest of the images due to cloud amount control, and this part of the images was not further processed in this study. We used three types of remote sensing sources, Landsat 8 OLI, Sentinel-2A, and their combined remote sensing factors, and then added topographic and climatic factors to model by RF. Finally, we compared them to pick the optimal model. The optimal model was ultimately used to predict the AGCS in the Pinus densata forest. The primary objective of this study was to explore a more appropriate method to estimate the AGCS in the Pinus densata forest in Shangri-La by remote sensing. Additionally, these efforts aim to fill the research gaps concerning the AGCS of the Pinus densata forest in Shangri-La.

2. Materials and Methods

2.1. Study Area

Shangri-La is located in the northwestern corner of Yunnan Province and is rich in biological resources. The area is situated in a highland area with a high terrain and an average elevation of approximately 3459 m. The terrain of Shangri-La is complex and varied, encompassing mountains, canyons, plateaus, grasslands, and lakes. Its administrative boundaries are bisected by the Jinsha River’s waterways in Diqing, Yunnan Province. The Jinsha River flows through the Shangri-La metropolitan area from the northwest to the southeast, while the Lancang and Nujiang Rivers flow through the area from the north and east, respectively. Shangri-La exhibits a subtropical monsoon climate alternating with a mountain monsoon climate. The region’s mean annual temperature is about 5.4 °C, and the annual precipitation is about 617 mm. The main tree species in the area include Pinus densata, Picea yunnanensis, Picea asperata, etc. [3]. In Shangri-La, Pinus densata participates in the maintenance of species diversity, regulates the climate, participates in the carbon cycle, and, as one of the dominant tree species in the region, accounting for 22.71% of the Shangri-La tree woodland area [50], has a huge potential for carbon sequestration. Pinus densata plays a significant role in the forestry industry of Shangri-La. The geographic location of the study area, the distribution of Pinus densata and the distribution of sample plots are shown in Figure 1.

2.2. Technical Route

In this research, we worked along the technical route shown in Figure 2, which mainly includes data processing; factor extraction, selection, and combination; modeling; and AGCS estimation.

2.3. Data Source and Processing

2.3.1. Remote Sensing Data

The remote sensing images we used are Landsat 8 OLI (30 m) and Sentinel-2A (10 m) in 2021, and they were downloaded from Geospatial Data Cloud (https://www.gscloud.cn/, accessed on 20 August 2023) and Copernicus Data Space Ecosystem (https://dataspace.copernicus.eu/, accessed on 20 August 2023), respectively. In order to ensure image quality, all remote sensing images we used had less than 1% cloud amount, so we did not use the cloud-masking operation. The details of each image data are shown in

Table 1. Information of the remote sensing images.

Image Type	Image ID	Cloud Amount/(%)	Acquisition Date
Landsat 8 OLI	LC81310412021082LGN00	0.88	23 March 2021
	LC81320412021313LGN00	0.58	9 November 2021
	LC81320402021313LGN00	0.89	9 November 2021
Sentinel-2A	S2A_MSIL2A_20211108T035951_N0500_R004_T47RNM_20230103T171747.SAFE	0.04	8 November 2021
	S2A_MSIL2A_20211118T040041_N0500_R004_T47RNK_20230101T032435.SAFE	0.05	18 November 2021
	S2A_MSIL2A_20211118T040041_N0500_R004_T47RNL_20230101T032435.SAFE	0.00	18 November 2021
	S2A_MSIL2A_20211118T040041_N0500_R004_T47RPK_20230101T032435.SAFE	0.17	18 November 2021
	S2A_MSIL2A_20211118T040041_N0500_R004_T47RPL_20230101T032435.SAFE	0.01	18 November 2021
	S2A_MSIL2A_20211218T040201_N0500_R004_T47RPM_20221225T151243.SAFE	0.00	18 December 2021

. All remote sensing images were used to extract the values on the image elements corresponding to the sample plots.

Landsat 8 OLI images were used after radiometric calibration and the FLAASH atmospheric correction. The Sentinel-2 images are level 2A, they were resampled to 10 m in SNAP 9.0.0, and mosaiced and cropped to use in ENVI 5.6. Finally, the resolution of Landsat 8 OLI images was downscaled to 10 m × 10 m by resampling in ENVI 5.6 [51]. Resampling of Landsat data may create some limitations, such as changes in spectral information in some areas, or result in changes or loss of some spatial information. We obtained three types of data sources.

• Source I: the factors from Landsat 8 OLI.
• Source II: the factors from sentinel-2A.
• Source III: the factors from the combination of Landsat 8 OLI and Sentinel-2A.

2.3.2. Sample Plots

The sample plots data were from 60 field plots surveyed in 2019 and 2021 [49]. All sample plots were 10 m × 10 m. All 60 sample plots, 20 from 2019 and 40 from 2021, were randomly distributed within the study area. The average AGB was calculated using the allometric growth equation of Pinus densata [52], which is also based on the measured average diameter at breast height (DBH) and tree height. Finally, the AGB density of each sample plot was obtained from the mean AGB and the area of the sample plots.

$$W=0.073·D^{1.793}·H^{0.880}$$

(1)

where W is the AGB per tree (kg), D is the diameter at breast height (cm), and H is the height (m).

The biomass density in the sample plots was multiplied by the Pinus densata carbon content coefficient to obtain the AGCS density. According to the Guidelines for Measuring Carbon Stock in Forest Ecosystems [53], the recommended carbon content coefficients apply to all major dominant tree species in China. The coefficient for Pinus densata is 0.501 [37]. Finally, the AGCS density of the sample plots was obtained by the carbon stock conversion formula [37].

$$C=W·CF$$

(2)

where C is the carbon stock (kg), W is the biomass (kg), and CF is the carbon content coefficient.

With Equations (1) and (2), we obtained the measured AGCS of Pinus densata in each sample plot as one of the key data for modeling.

2.3.3. Topography Data

According to the DEM data of Shangri-La, we extracted the topographic factors by ArcGIS, including elevation, slope, and aspect. The DEM data we used are from the ASTER GDEM data product. They were downloaded from Geospatial Data Cloud (https://www.gscloud.cn/, accessed on 20 August 2023). The data have a horizontal accuracy of 30 m and a vertical accuracy of 20 m. The DEM data were converted from 30 m to 10 m resolution by resampling. Thus, we obtained topographic factors with a resolution of 10 m.

2.3.4. Climate Data

Climatic factors were obtained from the National Tibetan Plateau Science Data Center, including mean annual temperature [54,55,56,57,58], annual precipitation [54,56,57,58,59], annual potential evapotranspiration, and monthly mean potential evapotranspiration [54,55,56,57,60,61] in 2021. The spatial resolution of this data is approximately 1 km. This dataset is generated by downscaling in China through the delta spatial downscaling scheme, which is based on the global 0.5° climate dataset released by CRU and the global high-resolution climate dataset released by WorldClim.

We used China’s monthly mean temperature data (0.1 °C) and monthly precipitation data (0.1 mm) [59] to calculate the mean annual temperature (1 °C) and annual precipitation (1 mm) by ArcGIS.

Similarly, we used the China’s potential evapotranspiration (0.1 mm) data [61] to calculate the monthly potential evapotranspiration (1 mm). The spatial resolution of this data is the same as that of temperature and precipitation data. It is obtained by using the Hargreaves potential evapotranspiration calculation formula [57,62,63] based on China’s 1 km monthly mean temperature, minimum temperature, and maximum temperature dataset [55].

2.4. Factors Extraction and Combination

2.4.1. Factors Extraction

To achieve data harmonization, we converted the Landsat 8 OLI images, DEM data, and climatic data from 1 km and 30 m to 10 m resolution by resampling. All the remote sensing variables in this work contain Landsat 8 OLI and Sentinel-2A texture factors, vegetation index factors, and band factors. Each texture feature can have some importance for AGCS modeling and estimation, but based on the results of factor selection, two kinds of texture features, CC and SM, are more important for AGCS modeling and estimation. The vegetation indices of Landsat 8 OLI were calculated by band math in ArcGIS on each sample plot. The vegetation indices calculation formula of Sentinel-2A was obtained from the SNAP. All remote sensing factors extraction for the two kinds of remote sensing images was performed in ArcGIS. Similarly, topographic and climatic factors were extracted by ArcGIS on each sample plot.

2.4.2. Factor Selection and Combination

The texture factors included eight types ranging from 3 to 11 window sizes in each band for both Sentinel-2A and Landsat 8 OLI images. The eight types of texture factors included homogeneity (HO), dissimilarity (DI), mean (ME), angular second-order moments (SM), entropy (EN), correlation (CC), variance (VA), and contrast (CO).

The band factors of Landsat 8 OLI images include B1~B7, B53, B64, B65, B67, B74, B547, and B4/Albedo. The vegetation indices factors include NDVI, TNDVI, RVI, SAVI, TSAVI, MASAVI, MSAVI2, GEMI, IPVI, and EVI.

The band factors of Sentinel-2A images include B1~B9, B11, and B12. The vegetation indices factors include NDVI, TNDVI, RVI, SAVI, TSAVI, MASAVI, MSAVI2, GEMI, IPVI, EVI, IRECI, MCARI, MTCI, REIP, NDI45, and PSSRa. The information of all extracted factors is presented in

Table 2. Information of all factors.

Types	Factors	Source
Topographic factors and band factors	Elevation, slope, aspect B1~B7/B1~B9, B11, B12 B53, B64, B65, B67, B74, B547, B4/Albedo	DEM, Landsat 8 OLI and sentinel-2A
Texture factors	(HO)homogeneity, (DI)dissimilarity, (ME)mean, (SM)angular second order moments, (EN)entropy, (CC)correlation, (VA)variance, (CO)contrast	Landsat 8 OLI and sentinel-2A
Vegetation indices factors	NDVI, TNDVI, RVI, SAVI, TSAVI, MASAVI, MSAVI2, GEMI, IPVI, EVI, IRECI, MCARI, MTCI, REIP, NDI45, PSSRa.	Landsat 8 OLI and sentinel-2A
Climatic factors	Mean annual temperature, annual precipitation, annual potential evapotranspiration, monthly mean potential evapotranspiration	National Tibetan Plateau Science Data Center

.

The results of factor selection are presented in

Table 3. Information of selected factors.

Data Type	Source of Data	Selected Factors
Source I	Landsat 8 OLI	LrR11B6CC, LrR11B5CC, LrR11B7CC, LrR11B6SM, LrR11B7SM
Source II	Sentinel-2A	SR5B8ASM, PSSRa, SR11B5SM, SR7B6CC, SR5B6CC, SR9B5SM, SR11B8ACC, SR5B1CC
Source III	Landsat 8 OLI and Sentinel-2A	LrR11B5CC, LrR11B6CC, LrR11B5SM, LrR11B7CC, LrR7B6CC, LrR9B6CC, SR5B6CC, SR5B8ASM, SR7B6CC, SR11B5SM, PSSRa

In the selected results, the expression is “S/LrRXBYZZ”. S/Lr is Sentinel-2A or resampled Landsat 8 OLI, RX is window size, BY is a certain band, ZZ is a certain texture feature’s abbreviation, and PSSRa is the chlorophyll index.

. The factors that have a significant correlation with AGCS density in the sample plots were selected by Pearson correlation analysis. Subsequently, we used RF to further select the factors with a cumulative feature importance contribution of 80% and above. For comparison, we performed the same selection of remote sensing factors of Landsat 8 OLI and sentinel-2A separately. VIF analysis is one of the commonly used multicollinearity detection methods [64,65,66], and the multicollinearity method for detecting the influence of unbiased predictive significance in RF modeling is average Gini impurity reduction (AGIR) comparison [67].

Another key data used for modeling was obtained by extracting factor values from selected remote sensing factors, topographic and climatic factors, through sample plots.

We added topographic and climatic factors in source I, II, III to model the following:

(1)

Combination of topographic and remote sensing factors, elevation and source I, II, III; slope and source I, II, III; aspect and source I, II, III.

(2)

Combination of climatic and remote sensing factors, annual precipitation and source I, II, III; mean annual temperature and source I, II, III; annual potential evapotranspiration and source I, II, III; monthly mean potential evapotranspiration and source I, II, III.

With all the combinations of data and factor additions, we obtained twenty-four sets of data for modeling:

Model 1: established by source I

Model 2: established by source II

Model 3: established by source III

Model 4: established by elevation and source I

Model 5: established by slope and source I

Model 6: established by aspect and source I

Model 7: established by elevation and source II

Model 8: established by slope and source II

Model 9: established by aspect and source II

Model 10: established by elevation and source III

Model 11: established by slope and source III

Model 12: established by aspect and source III

Model 13: established by annual precipitation and source I

Model 14: established by mean annual temperature and source I

Model 15: established by annual potential evapotranspiration and source I

Model 16: established by monthly mean potential evapotranspiration and source I

Model 17: established by annual precipitation and source II

Model 18: established by mean annual temperature and source II

Model 19: established by annual potential evapotranspiration and source II

Model 20: established by monthly mean potential evapotranspiration and source II

Model 21: established by annual precipitation and source III

Model 22: established by mean annual temperature and source III

Model 23: established by annual potential evapotranspiration and source III

Model 24: established by monthly mean potential evapotranspiration and source III

2.5. Model Establishment and Evaluation

In this research, we used RF to model and estimate the AGCS in the Pinus densata forest. This method is primarily utilized to solve classification and regression problems. It has the advantages of high accuracy, resistance to overfitting, and interpretability. The RF parameters were described in Bao et al. [47]. These parameters and their ranges in this work are: n_estimators: 20–200, max_depth: 6–9, min_samples_leaf: 2, and min_samples_split: 2. These parameters were chosen based on multiple cross-validation tests. The dataset for this study is small, so we used the parameters described above. Using too large a range of parameters may result in overfitting.

We used RF to randomly select 80% of the sample data for model training. The remaining 20% of the sample data was used to validate the trained model. The model evaluation indicators include the coefficient of determination (R²), the root mean square error (RMSE), the relative root mean square error (rRMSE), and accuracy (P). The R² is used to reflect the overall model fitting effect, the RMSE is used to quantify the degree of error, the rRMSE is used to indicate the relative deviation of the estimated value from the measured value, and P is used to reflect the average prediction accuracy of the model [68]. The formulas are as follows:

$$R^2={\displaystyle \frac{{\sum }_{i=1}^n{\left(\widehat{y_i}-\overline{y}\right)}^2}{{\sum }_{i=1}^n{\left(y_i-\overline{y}\right)}^2}}$$

(3)

$$RMSE=\sqrt{{\displaystyle \frac{{\sum }_{i=1}^n{\left(y_i-\overline{y}\right)}^2}{n}}}$$

(4)

$$rRMSE={\displaystyle \frac{RMSE}{\overline{y}}}$$

(5)

$$P={\displaystyle \frac{1}{n}}{\sum }_{i=1}^n\left(1-\left|{\displaystyle \frac{y_i-{\widehat{y}}_i}{{\widehat{y}}_i}}\right|\right)$$

(6)

where $y_i$ is the true value, $\widehat{y_i}$ is the model regression value, $\overline{y}$ is the mean value, and $n$ is the plot number.

We compared all the models and picked an optimal model. Finally, we used the data used in this model to perform AGCS estimation by RF as a way to develop a method for AGCS estimation based on remote sensing data and topographic or climatic data.

3. Results

3.1. Modeled by Remote Sensing Factors

Source I, source II, and source III were modeled before adding topographic or climatic factors. The results are shown in

Table 4. Modeling results from different data sources.

Data Source	Model	R2	RMSE/(t·ha−1)	rRMSE/(%)	P/(%)
Source I	Model 1	0.85	11.38	23.46	78.71
Source II	Model 2	0.82	12.41	24.21	79.74
Source III	Model 3	0.87	10.81	23.19	79.71

.

Among the three models, the model established by source III has the best R², RMSE, and rRMSE. Its indicators are better than those of the models constructed by source I and source II. The accuracy of the model based on source III is also very close to that of the model based on source II. From the modeling results in Table 4, it can be seen that the source I-based model has the advantage of better R² and RMSE, while the source II-based model has the advantage of better rRMSE and P. The model based on source III has more balanced model indicators, which may be due to the fact that the high resolution of the Sentinel-2A data provides more detail and information. The scatterplot of each model is shown in Figure 3.

From a practical application point of view, the coefficient “b” can be interpreted as a baseline adjustment term. In biomass estimation research, due to the complexity and uncertainty of the measurement process, it may not be possible to fully and accurately describe the real situation of biomass only by relying on the linear relationship between the measured value and the predicted value. At this time, the coefficient “b” can fine-tune this linear relationship to make the predicted value closer to the actual biomass.

3.2. Modeled Adding Topographic Factors

The topographic factors were added to source I, source II, and source III for modeling. The results are shown in

Table 5. Modeling results adding topographic factors.

Data Source	Added Factors	Model	R2	RMSE/(t·ha−1)	rRMSE/(%)	P/(%)
Source I	Elevation	Model 4	0.88	10.37	21.90	80.67
	Slope	Model 5	0.88	10.24	21.57	80.92
	Aspect	Model 6	0.88	10.38	22.33	81.14
Source II	Elevation	Model 7	0.85	11.51	21.49	82.70
	Slope	Model 8	0.86	11.22	21.69	82.10
	Aspect	Model 9	0.85	11.37	20.56	82.80
Source III	Elevation	Model 10	0.89	10.01	21.47	82.17
	Slope	Model 11	0.88	10.16	21.34	81.33
	Aspect	Model 12	0.89	9.92	21.95	81.18

. A comparison of the modeling results after adding topographic factors is shown in Figure 4.

As shown in Figure 5, after adding the topographic factors, various indicators of the AGCS estimation model have different improved percentages. The comparative results show that the R² of the model based on source III is the highest, the R² of the model based on source I is the second highest, and the R² of the model based on source II is the lowest group after adding the same topographic factor. Regarding prediction accuracy, the model based on source II is the highest when adding the same topographic factor; after adding aspect, the prediction accuracy of model 9 reaches 82.80%. The model based on source III is the second highest; after adding elevation, the prediction accuracy of model 10 reaches 82.17%. The prediction accuracy of the model based on source I is the lowest relatively; after adding aspect, the prediction accuracy of model 6 reaches 81.14%. Furthermore, the model indicators improved their percentage when adding the topographic factor in the order of model source II, model source I, and model source III.

3.3. Modeled Adding Climatic Factors

As previously stated, the climatic factors were added to source I, source II, and source III to model. After that, the modeling results were subjected to a comparative analysis. The model indicators are presented in

Table 6. Modeling results after adding climatic factors.

Data Source	Added Factors	Model	R2	RMSE/(t·ha−1)	rRMSE/(%)	P/(%)
Source I	AP	Model 13	0.88	10.34	21.51	81.54
	MAT	Model 14	0.89	9.81	20.24	82.59
	APET	Model 15	0.89	10.02	20.16	82.83
	MMPET	Model 16	0.88	10.40	21.26	82.70
Source II	AP	Model 17	0.88	10.77	19.60	84.42
	MAT	Model 18	0.85	11.40	20.66	83.11
	APET	Model 19	0.86	11.18	20.69	82.91
	MMPET	Model 20	0.85	11.41	21.83	83.00
Source III	AP	Model 21	0.90	9.53	20.59	83.00
	MAT	Model 22	0.89	9.95	21.39	81.82
	APET	Model 23	0.88	10.07	21.54	82.03
	MMPET	Model 24	0.89	9.88	22.02	81.85

. Figure 6 illustrates the comparative analysis of the modeling improvement after adding the climatic factors.

As illustrated in Figure 7, the modeling results after adding climatic factors to the three data sources are as follows: all the model indicators exhibited better values in comparison to those of the model constructed before adding topographic and climatic factors. Compared with the modeling results when adding topographic factors, the R² and prediction accuracy of models 13–24, which included the climatic factors, has a further improvement. Among the three data sources, the R² of the model based on source III is once again the highest. The R² of the model based on source I is the second highest, but it is very close to the R² of the model based on source III. The R² of the model based on source II is the lowest relatively. Among the three data sources, the prediction accuracy of the model based on source II including climate factors remains the highest; after adding annual precipitation, the prediction accuracy of model 17 reaches 84.42%. The prediction accuracy of the model based on source III is the second highest, after adding annual precipitation, the prediction accuracy of model 21 reaches 83.00%. Among the three data sources, the prediction accuracy of the model based on source II is the lowest relatively. The order of improved percentage of model indicators is source II, source II, and source III, in descending order. A comparison of the model established after adding factors to source III indicates that the most effective climatic factor is the annual precipitation, and the most effective topographic factor is the elevation.

3.4. AGCS Mapping

We compared the AGCS estimation models of Pinus densata forest, which were established by adding topographic or climatic factors to the three kinds of data sources. Then, the optimal model obtained by adding the annual precipitation factor to source III was used for AGCS estimation and mapping, and it has been shown that Pinus densata is more sensitive to moisture effectiveness (precipitation and relative humidity) in the early part of the growing season [62]. Finally, the estimated AGCS value for the Pinus densata forest in Shangri-La was 9,737,500 t in 2021. The distribution of AGCS is illustrated in Figure 8. Some scholars have also made remote sensing estimates of Pinus densata AGB and AGCS in Shangri-La using different methods and data. Compared to the following studies, we have not only added some climatic factors, but also Sentinel-2A data. The carbon stock of the Pinus densata arbor layer in 2008, estimated by Yue [50], was 8,640,900 t based on Landsat TM image data. Wang et al. [45] used remote sensing information modeling to estimate the AGB in 2009, and the result was 20,000,000 t. Liao et al. [29] estimated the AGB during the 1987–2017 period by Landsat time series images and RF algorithm, and the result was 8,496,300~9,157,800 t. Teng et al. [8] combined the AHTC filtering algorithm and RF algorithm to estimate the AGB during 1987–2017 similarly; the results ranged from 11,545,300 t to 16,542,000 t. Sun [31] estimated the AGB by Landsat 8 OLI data, and the result was 11,719,600 t. Xie [69] used Landsat 8 OLI images combined with the K-NN algorithm to estimate the AGB in 2015, and the result was 12,100,000 t. Chen et al. [49] based their analyses on RF and combined Sentinel-1 and Sentinel-2A data to estimate the AGB in 2021, and the result was 17,210,000 t. In addition to the results of the AGCS study conducted by Yue [50], we calculated the AGB of the studies mentioned above to AGCS values using the Pinus densata carbon content coefficient of 0.501. These values are presented in

Table 7. Biomass and carbon stock values estimated by relevant studies.

Data Year	AGB Value/Million Tons	AGCS Value/Million Tons	Source	Methodology
2008	16.67	8.64	Yue [50]	Modeled and estimated AGB by SVM
2009	20.00	10.02	Wang et al. [45]	Combining multiple factors to estimate AGB by remote sensing information model
2015	12.10	6.06	Xie [69]	Adding topographic factors and estimating AGB by optimized k-NN
2016	11.72	5.87	Sun [52]	Establishing the biomass model of sample trees and then estimating overall AGB
1987–2017	8.50~9.16	4.26~4.59	Liao et al. [29]	Adding topographic factors and establishing dynamic model to estimate AGB by RF
1987–2017	11.55~16.54	5.78~8.29	Teng et al. [8]	Estimating AGB after improving image quality by filtering algorithms
2021	17.21	8.62	Chen et al. [49]	Estimating AGB by Sentinel-1/2 data and RF while considering the seasonal effects
2021	\	9.74	This study	Estimating AGCS by RF based on two kinds of remote sensing data and climatic factors

. Compared with the research results of Chen et al. [49] and Xie et al. [69], the AGCS values of the three studies are relatively close, and our estimation result is not the highest one. Therefore, our results are reliable.

According to Table 7, the AGCS values estimated in this study are in between the above results. The carbon stock values of the related studies listed in the table fluctuate up and down. This phenomenon may be attributed to some reasons. These reasons include the types of remote sensing images used in each study, the sources and errors of the sample plots, the methods of calculating carbon stock, the differences in modeling methods, the characteristics of the algorithms utilized for modeling, and the differences in the Pinus densata forest area between periods, which collectively result in fluctuating and imprecise values in the estimation. In addition, the sample plot data used in this study were surveyed in 2019 and 2021, respectively, and tree growth during this period may also affect AGCS estimation. However, the estimated annual growth of Pinus densata in the region is very small, and the impacts would be minimal [49].

We estimated the AGCS value of each distribution point based on the Pinus densata distribution data of the corresponding year, the optimal model data, and the RF algorithm to get the final AGCS distribution map in ArcGIS 10.8.1.

4. Discussion

4.1. Application of Remote Sensing Data Combination in Forest AGCS/AGB Estimation

Some studies of forest AGB or AGCS estimation have selected only Landsat as the remote sensing data [13,24,25,26,27,28,29,30,31,32,37]. Single Landsat data may bring some insufficiency. According to this research, although the models based on Landsat have a higher R², their prediction accuracy is relatively low [8,29,37,47]. Other studies used only Sentinel-2 data [15,16,17,18,19,20,23,48]. Although Sentinel-2 data can provide rich spectral information, the information provided by a single remote sensing data is limited. In addition, the model fitting effect based on Sentinel-2 data is poor relative to that based on Landsat [69,70,71]. In this study, we also found that the AGCS estimation models constructed by Landsat 8 OLI data generally have better fitting effects, while the AGCS estimation models established by Sentinel-2A data have better prediction accuracy. Combining the two kinds of remote sensing image data, we can obtain more spectral information, texture factors, and vegetation indices. Therefore, we obtained an AGCS estimation model with excellent fitting effect and prediction accuracy. In the other studies [23,49,50,72,73,74,75,76,77,78], the authors combined different remote sensing data to study the AGB and stock volume. These data included Landsat 8 OLI, Sentinel-1, Sentinel-2, LiDAR, and others. It can be seen from these studies that more remote sensing data can bring more information. Moreover, there are AGB multi-source remote sensing estimation studies that combined UAV Lidar data and UAV multispectral data [79]. There is also a biomass estimation study based on multi-source remote sensing data combined with deep learning [80] and AGB multi-source remote sensing estimation studies using a combination of GF-2, Landsat 8, and Airborne Laser Scanning (ALS) data [81]. At present, the application of multi-source remote sensing is increasing. Compared to these studies, we lack LiDAR data or SAR data, but environmental factors also need to be taken into account in AGB estimation studies. Deep learning may also receive more attention in future AGB multi-source remote sensing estimation studies. However, the reason why we only used two types of remote sensing data is that this study aims to explore the influence of adding factors on the estimation of Pinus densata AGCS. Therefore, we did not use more complex multi-source remote sensing data. In addition, differences in sensors or spectral ranges may also affect the comparability of AGB and AGCS estimation based on multi-source remote sensing data.

4.2. Advantages in Model Accuracy from Sentinel-2A

According to modeling results, the prediction accuracy of the models based on source II is always the best among the three kinds of data sources. It aligns with the findings of Zhou and Feng [74], which demonstrated that the model based on the Sentinel-2 images achieved better accuracy compared to the same model based on the Landsat 8 factors set. Puliti et al. [75] also mentioned that Landsat data had lower prediction accuracy than Sentinel-2 data. This may be because the Sentinel-2A product has a relatively large spatial coverage and high resolution, which is more advantageous to the estimation of the AGB and AGCS [78]. This is consistent with the result mentioned by Huang et al. [76]: the prediction accuracy is higher for Sentinel-2 data than Landsat 8 data. The models based on source III arrived at the second highest prediction accuracy, which may be related to the resolution of the remote sensing data itself and the size of sample plots. The Sentinel-2 images resolution supports an accuracy of up to 10 m, while the Landsat 8 OLI images have a resolution of up to 30 m. Combining two kinds of remote sensing factors may have a complementary effect on the model fitting effect and prediction accuracy, and the model based on source III has satisfactory results in both the fitting effect and prediction accuracy. This is consistent with the findings of Luo et al. [77], who demonstrated that the combination of Sentinel-2B and Landsat 8 OLI data yielded superior results. Regarding the modeling results, the model based on source III already had good results before the addition of topographic or climatic factors. Although these factors have a limited degree of optimization for the model, the improvement of the models based on source I and source II is much higher, while the improvement of the models based on source III was relatively poorer. Nevertheless, when adding annual precipitation, the model prediction accuracy based on source III reached 83.00%. It was higher than some previous studies. For instance, the prediction accuracy of AGB was 78.77% ± 2.39% with the addition of topographic factors by Liao et al. [29]. The prediction accuracy of carbon stock was 81.00% with the addition of land surface temperatures by Yin et al. [37].

4.3. Changes in Model Indicators and the Importance of Climatic Factors

From Table 4, Table 5 and Table 6, it can be seen that the changes in indicators between models based on the same data source are not obvious. It is due to the fact that the model indicators represent the relationship between the measured values and the predicted values [68], and the predicted values usually do not have a large gap with the measured values, which results in the indicators of the models not changing much as well. If the indicators of each model change too much, it means that the predicted values of the model change a lot, and the model is unstable. In addition, the multicollinearity between remote sensing factors may also cause insignificant differences in indicators between models.

In summary, it is difficult to get a significant change in model metrics based on the same data sources, but the improvement in model indicators from the models based on source I and II to the models based on source III when adding annual precipitation is more significant. We also demonstrated that, in some regions, climate factors have a positive role in AGB or AGCS estimation. In future studies, local climate factors can be considered to get higher accuracy of AGB and AGCS estimation.

Topographic factors are often used in studies estimating forests AGB and AGCS [24,25,26,27,28,29,30,31,32]. However, considering only topographic factors is not comprehensive enough. Without major natural disasters, the topography of these areas remains unchanged for long periods [82]. Therefore, the effect of topography on the AGCS is relatively constant. However, the light and rainfall that the trees receive are affected by topography [80]. It also affects the climatic conditions where the trees are located.

At present, climate change is occurring on Earth. Therefore, we added some climatic factors to the estimation of Pinus densata AGCS and obtained satisfactory results. As is known to all, forest ecosystems have a certain degree of influence on the global climate. However, climate does not unilaterally affect the forest ecosystem. Instead, the two influence each other. A study [83] has shown that the carbon density of Larix principis-rupprechtii is affected by temperature and precipitation. Moreover, mean annual temperature and annual precipitation are essential factors affecting carbon stock in forest ecosystems [42]. Another study [84] has shown that mean annual temperature and annual precipitation have different degrees of influence on the growth and carbon sequestration capacity of five primary planted forests (Larix. spp., Pinus massoniana, Cunninghamia lanceolata, Populus spp., Eucalyptus spp.) in China. Regarding carbon stock, the carbon stock changes of the Atlantic forest biome are more sensitive to mean annual temperature and annual precipitation [85]. This is why we think that under certain conditions, the topographic factors are not so important. We believe that if climate change is evident in a region, it is more important to consider climatic factors when estimating the AGCS in that region.

From Equation (1) [52], we can know the direct relationship between the AGB of Pinus densata and DBH and tree height, which are important reference indexes for tree growth. The effects of temperature and precipitation on the growth of Pinus densata have been studied by some scholars. One of the studies stated that climate warming and drying and uneven seasonal precipitation lead to growth limitation of Pinus densata [86]. There are also studies indicating that the relationship between Pinus densata growth and climate is variable [62] and that Pinus densata in different places respond differently to climatic factors but are all affected by temperature and moisture [87]. There are also studies indicating that Pinus densata growth is mainly limited by water availability early in the growing season [88].

It is easily understood that climatic factors are significant for AGCS estimation and should be considered in relevant studies. We also hope that more accurate climate data will contribute to the future estimation of forest carbon stock.

5. Conclusions

This study is based on Landsat 8 OLI and Sentinel-2A images. We added topographic or climatic factors to the selected remote sensing factors and modeled by RF, and then compared model prediction accuracy to select the optimal model. The optimal model is used to estimate the AGCS. The following conclusions were obtained: (1) texture factors are more important for AGCS modeling and estimation than vegetation indices (2) adding topographic or climatic factors can improve the prediction accuracy of the AGCS; (3) adding climatic factors improves the model accuracy more than adding topographic factors; and (4) of all the climatic factors, the addition of the annual precipitation factor provided the greatest improvement in model prediction accuracy. We developed a more accurate AGCS estimation method based on Landsat 8 OLI data, Sentinel-2A data, and annual precipitation data. This research’s results can provide a reference for forest AGCS estimation based on Landsat 8 OLI and Sentinel-2A data. Higher-resolution climate data will be needed in future AGCS estimation studies.