The Carbon Neutral Potential of Forests in the Yangtze River Economic Belt of China

Abstract

Prediction of forest carbon sink in the future is important for understanding mechanisms concerning the increase in carbon sinks and emission reduction, and for realizing the climate goals of the Paris Agreement and global carbon neutrality. Based on stand volume data of permanent monitoring plots of the successive national forest inventories from 2004 to 2018, and combined with multiple variables, such as climatic factors, soil properties, stand attributes, and topographic features, the random forest algorithm was used to predict the stand volume growth-loss and then calculated the forest biomass and its carbon sink potential between 2015 to 2060 in the Yangtze River Economic Belt of China. From 2015 to 2060, the predicted forest biomass carbon storage and density increased from 3053.27 to 6721.61 Tg C and from 33.75 to 66.12 Mg C hm⁻², respectively. The predicted forest biomass carbon sink decreased from 90.58 to 73.98 Tg C yr⁻¹, and the average forest biomass carbon storage and sink were ranked in descending order: Yunnan, Sichuan, Jiangxi, Hunan, Guizhou, Hubei, Zhejiang, Chongqing, Anhui, Jiangsu, and Shanghai. The forest biomass carbon storage in the Yangtze River Economic Belt will increase by 3.67 Pg C from 2015 to 2060. The proportion of forest C sinks on the regional scale to C emissions on the national scale will increase from 2.9% in 2021–2030 to 4.3% in 2041–2050. These results indicate higher forest carbon sequestration efficiency in the Yangtze River Economic Belt in the future. Our results also suggest that improved forest management in the upper and middle reaches of the Yangtze River will help to enhance forest carbon sink in the future.

1. Introduction

Since the Industrial Revolution, anthropogenic greenhouse gas emissions (e.g., CO₂, CH₄, and N₂O) from the consumption of fossil fuels have accelerated global warming, which poses huge challenges to terrestrial ecosystems and human well-being [1]. Carbon (C) sequestration in terrestrial ecosystems plays a vital role in mitigating global warming [2,3,4]. Over the past 50 years, terrestrial ecosystems have been responsible for removing approximately one-third of anthropogenic CO₂ emissions [5]. As the largest C pools in the terrestrial ecosystem, above- and below-ground components of the forest ecosystems comprise more than 80% and 40% of global terrestrial C storage, respectively [6]. Therefore, the prediction of forest C sinks in the future is important for understanding the mechanisms of the increase in C sinks and emission reduction, and for realizing the climate goals of the Paris Agreement and global C neutrality.

China recently announced a long-term climate mitigation plan to peak emissions before 2030 and be C neutral by 2060. Many innovative technologies, such as renewable energy production, food system transformation, waste valorization, C sink conservation, and C-negative manufacturing, can offer solutions for achieving C neutrality [7]. Compared to technical potentials, cost-effective (available up to USD 100/t CO₂ eq) estimates represent a more realistic and actionable target for policy and are approximately 50% of forests and other ecosystems [4]. Therefore, C removals using natural C sinks in terrestrial ecosystems are vital [4]. China’s forested area (especially plantations) is increasing [8], and China is one of the top countries with the highest total cost-effective mitigation potential for global warming from land-based measures [4]. Accordingly, several studies for estimating forest C storage in China have been carried out to evaluate forest C neutrality potential [9,10,11,12]. However, the forest C sink prediction in previous studies was often analyzed based on the stand attributes (e.g., age and biomass increment) and several climatic factors (temperature and precipitation) [9,12]. Other studies also found that forest C sink was affected by the edaphic properties (mainly nitrogen and phosphorus) and the topographic features (e.g., altitude and aspect) [13].

Moreover, afforestation is considered a cost-effective and readily available climate change mitigation option [14]. China has significant experience with large-scale afforestation and reforestation programs [15,16], which resulted in a 25% net increase in global canopy area on 6.6% of the global vegetated area between 2000 and 2017 [17]. But a recent estimate reflects a previously underestimated land C sink over specific regions (e.g., southwest China) throughout the year [18]. Therefore, taking possible climatic factors, soil properties, stand attributes, and topographic features into account, the way to accurately estimate forest C storage in China across scales and regions needs to be better understood.

The Yangtze River Economic Belt (YREB), including 11 provinces and municipalities (

Table A1. Variation in the forest area of 11 provinces and municipalities in the Yangtze River Economic Belt based on the National Forest Inventory of China (10⁶ hm²).

Province	National Forest Inventory
	6th (1999–2003)	7th (2004–2008)	8th (2009–2013)	9th (2014–2018)
Shanghai	0.02	0.06	0.07	0.09
Jiangsu	0.77	1.08	1.62	1.56
Zhejiang	5.54	5.84	6.01	6.05
Anhui	3.32	3.60	3.80	3.96
Jiangxi	9.31	9.74	10.02	10.21
Hubei	4.98	5.79	7.14	7.36
Hunan	8.61	9.48	10.12	10.53
Chongqing	1.83	2.87	3.16	3.55
Sichuan	14.64	16.60	17.04	18.40
Guizhou	4.20	5.57	6.53	7.71
Yunnan	15.60	18.18	19.14	21.06
Total	68.83	78.80	84.66	90.48

), is the most important economic development and ecological protection region in China, indicating the coexistence of high CO₂ emissions and high forest C sinks [19,20]. Although the land area only accounts for 21.9% of China’s total land area, the forest area accounts for 41.0% of the total forest area [8]. As a result, forest ecosystems in this region have a high potential for C sink. Owing to the high heterogeneity of topographic features and soil fertility (Figure 1 and

Table A2. Candidate variables for stand volume growth-loss modeling.

Type	Variable	Mean	SD	Min.	Max.	Type	Variable	Mean	SD	Min.	Max.
Stand characteristics	Stand age (yr)	31.10	34.86	0.00	300.00	Climatic factors	Bio1 (°C)	16.33	2.35	3.50	20.32
	Volume gross growth (m3·ha·yr−1)	4.97	3.18	0.00	16.00		Bio2 (°C)	9.70	2.28	4.77	19.18
	mortality (m3·hm−2·yr−1)	1.88	2.60	0.00	9.68		Bio3 (°C)	31.20	2.62	21.59	34.38
	Dominant tree species	-	-	-	-		Bio4 (°C)	1.79	3.26	−15.22	9.57
Topography	Altitude (m)	1082.28	1063.87	0.00	5000.00		Bio5 (°C)	29.21	3.55	14.82	35.69
	Slope direction	-	-	-	-		Bio6(°C)	23.54	2.57	11.61	26.42
	Slope position	-	-	-	-		Bio7 (°C)	12.21	2.62	−1.48	17.88
	Gradient (°)	24.54	12.13	0.00	80.00		Bio8 (°C)	24.91	3.10	12.06	28.24
Soil variables	Soil thickness (cm)	59.82	23.04	1.00	300.00		Bio9 (°C)	6.82	2.93	−5.66	14.96
	AN (mg·kg−1)	133.63	71.33	18.48	640.75		Bio10 (°C)	1223.88	311.47	556.41	1919.96
	AP (mg·kg−1)	5.28	3.00	0.93	44.29		Bio11 (mm)	231.13	46.40	109.78	376.93
	AK (mg·kg−1)	122.99	57.35	22.28	377.47		Bio12 (mm)	27.16	19.79	0.93	64.91
	BD (g·cm−3)	1.22	0.15	0.47	1.45		Bio13 (mm)	538.85	117.33	256.26	814.22
	PH	5.83	1.01	4.30	9.32		Bio14 (mm)	152.45	50.33	35.98	247.20
	GRAV (%)	8.68	7.56	0.02	48.73		Bio15 (mm)	554.72	84.91	278.09	853.31
	SOM (g·kg−1)	3.86	2.58	0.39	33.26		Bio16 (mm)	112.58	84.00	1.29	300.99

Notes: (1) “-” represents categorical variables; SD, standard deviation; Min, minimum; Max, maximum. (2) Stand age, average age of dominant tree species in main layer of arbor forest; Volume gross growth, average annual volume gross growth of forest stand; Volume loss, average annual volume loss of forest stand; AN, alkali-hydrolysable N; AP, available P; AK, available K; BD, bulk density; PH, soil pH value in the topsoil; GRAV, gravel; SOM, soil organic matter; Bio1, annual mean temperature; Bio2, mean diurnal range (mean of monthly [max temp.-min temp.]); Bio3, max temperature of warmest month; Bio4, min temperature of coldest month; Bio5, temperature annual range; Bio6, mean temperature of wettest quarter; Bio7, mean temperature of driest quarter; Bio8, mean temperature of warmest quarter; Bio9, mean temperature of coldest quarter; Bio10, annual precipitation; Bio11, precipitation of wettest month; Bio12, precipitation of driest month; Bio13, precipitation of wettest quarter; Bio14, precipitation of driest quarter; Bio15, precipitation of warmest quarter; Bio16, precipitation of coldest quarter.

), however, the forest C storage estimates are uncertain [21]. Additionally, while previous studies have estimated the vegetation C storage in the YREB [22], they have not included the new afforestation forests with a high potential for climate change mitigation [14]. Therefore, more accurate estimates of forest C storage in the YREB are needed.

In this study, based on the national forestry inventory (NFI) data of volume and environmental factors (e.g., climate, soil, topography, and stand attributes), we performed the random forest (RF) algorithm to construct the growth-loss model of the dominant tree species (or groups) and then to predict forest biomass and its C storage from 2015 to 2060 in the YREB. The main objectives of this study were (1) to accurately estimate the future forest biomass C sink potential and its spatial-temporal dynamic, and (2) to explore the contribution of forest biomass C sink in this region to CO₂ emission reduction in China. We expect that this study will provide scientific guidance for accurately positioning the role of forests in the YREB in increasing sinks and reducing emissions, and scientifically formulating China’s “C peak and C neutrality” action path and target management.

2. Materials and Methods

2.1. Study Area

The YREB is a major national strategic development area (both economy and ecology) in China which spans the three major regions of China’s middle, east and west, within longitudes of 90°31′50″–121°53′23″ E and latitudes of 21°8′45″–34°56′47″ N (Figure 1). This region covers 11 provinces and municipalities (Table A1), and its total area is roughly 2.05 million km² [22]. There is diverse topography, soil type, land use, and land cover [22], and the forest is the largest vegetation type in the YREB. Forest area is increasing from 1999 to 2018 (Table A1). The forest area is 90.48 million hm^2, and its coverage is 44.38%, accounting for 41.04% of the total Chinese forest in the 9th NFI (Table A1). From east to west, elevation gradually rises with a maximum difference of 7000 m (Figure 1), leading to the Yangtze River flowing west to east. Moreover, except for the Plateau climate of the western region, most regions are typical of subtropical monsoon climate [21,22]. Mean annual temperature (MAT) is between 3.50 °C to 20.32 °C, and mean annual precipitation (MAP) is between 556.41 mm to 1919.96 mm. Soil types mainly include Haplic Acrisol soil and Haplic Luvisol soil [22].

2.2. Data Sources

The NFI database includes detailed information on forest vegetation and sampling plot [12]. The stand volume data were collected from the permanent monitoring plots during the 7th (2004–2008), 8th (2009–2013), and 9th (2014–2018) period of NFI in China. There were 21,433 permanent plots of arbor forests, 4646 permanent plots of shrub, and 1515 permanent plots of bamboo in the YREB. The plots are generally square with an area of 0.067 hm² and are monitored every 5 years. The site conditions (e.g., altitude, slope, aspect, slope position, soil thickness, land classification) and stand attributes (e.g., origin, dominant tree species, stand age, volume gross growth, mortality, and volume) of each plot were collected (Table A2). The stand volume is the total volume of all standing trees. The stand volume growth, in terms of volume increment per hectare per year, equals the stand volume difference between the two surveys. The descriptive characteristics of these indicators are listed in Table A2.

In this study, seven edaphic variables were chosen as candidate variables for model fitting, including alkali-hydrolyzable nitrogen (AN), available phosphorus (AP), available potassium (AK), bulk density (BD), pH value (PH), gravel content (GRAV), and soil organic matter (SOM) (Table A2), which were derived from the data of “Soil properties of land surface in China” [23] released by the platform of Tibet Plateau Science Data Center (http://tpdc.ac.cn). Based on 8979 soil profiles and the soil map of China (1:1,000,000), the database was established using the polygon linkage method with a spatial resolution of 30 × 30 arc-seconds. It has been widely used in related research [24,25]. In this study, we assumed that all edaphic variables would not change in the future.

Bioclimatic factors include historical climate data used for modeling and future climate scenario data used for prediction. The historical climate data comes from temperature and precipitation databases between 1995 and 2015 provided by China’s Meteorological Data Sharing System (http://data.cma.cn). The national meteorological factor raster layer data is from the Co-Kerry Golden interpolation (Kring). The future climate data comes from the Shared Socio-economic Pathways of 370 in the global climate models of MIROC-ES2L of the sixth Coupled Model Intercomparison Project (CMIP6) in the World Climate Database (http://www.worldclim.org) from 2021 to 2060. Based on the latitude and longitude information of each permanent plot, historical and future scenario climate data in these plots were obtained from these databases, including 16 climatic factors (Table A2).

2.3. Forecast of Forest Area in the Future

The Forest area of the YREB in the future includes two parts: the existing forest area and the new afforestation area. The existing forest area referred to the forest area of the 9th NFI from 2014 to 2018 (hereafter 2015). In this study, we assumed that the existing forest area would not be affected by anthropogenic disturbances such as forest conversion, deforestation, and logging in the future.

Combined with the data of suitable areas for forest growth in China [26], land cover data [27], and DEM data (Geospatial Data Cloud site, http://www.gscloud.cn), the distribution area of potential new afforestation areas across China was launched after deducting the existing forest land, cultivated land with a slope of less than 25°, wetlands, water bodies, tundra, impervious surfaces, and ice-snow. Based on this, the total area of potential new afforestation in the provinces and municipalities in the YREB was determined. China’s forestry planning stated that the “forest coverage rate will increase to 23.04% by 2020, 24.1% by 2025, and then stabilize above 26% by 2050”. Therefore, based on the changes in land types during the 6th to 9th NFI, the area of new plantations and natural regeneration forests in each province during this period is calculated. It is assumed that the area proportion of new plantations and natural regeneration forests in the future will be the same as those of the 6th to 9th NFI, so as to determine the growth of forest area in the provinces and municipalities of the YREB in the future (

Table A3. Afforestation and reforestation area of forest from 2016 to 2060 in YREB (10⁴ hm²).

Province	Forest Type	Periods
		2016–2020	2021–2030	2031–2040	2041–2050	2051–2060
Jiangsu	Arbor	0.2	0.35	0.33	0.22	0
Zhejiang	Arbor	6.27	11.35	10.38	6.9	0
	Bamboo	1.32	2.39	2.19	1.46	0
Anhui	Arbor	2.13	3.85	3.52	2.34	0
	Bamboo	0.27	0.49	0.45	0.29	0
Jiangxi	Arbor	5.77	10.44	9.54	6.35	0
	Bamboo	0.75	1.36	1.24	0.83	0
Hubei	Arbor	8.11	14.69	13.42	8.92	0
	Bamboo	0.24	0.43	0.39	0.26	0
Hunan	Arbor	13.05	23.62	21.6	14.36	0
	Bamboo	1.34	2.43	2.22	1.48	0
Sichuan	Arbor	56.02	101.43	92.7	61.67	0
	Shrub	0.01	0.02	0.02	0.01	0
	Bamboo	2.49	4.51	4.13	2.74	0
Guizhou	Arbor	33.19	60.09	54.92	36.54	0
	Bamboo	0.91	1.64	1.5	1	0
Yunnan	Arbor	62.21	112.62	102.94	68.48	0
	Bamboo	0.38	0.7	0.63	0.42	0
Chongqing	Arbor	5.94	10.76	9.83	6.54	0
	Bamboo	0.37	0.68	0.61	0.41	0

) until the upper limit of the total potential new forest area is reached.

2.4. Model Development, Validation, and Prediction

This study assumed that the biomass per unit area of bamboo and shrub forests remained stable; that is, the annual growth was zero. The change in total biomass of bamboo and shrub forests was only related to the area change. This study only used the random forest algorithm to simulate and predict the future growth of arbor forests.

2.4.1. Random Forest Algorithm

The primary modeling tool used in this study is the Random Forest (RF) model, an ensemble method based on decision trees that consider the ability to solve regression problems and classification problems [28]. The model constructs a series of base learners through resampling and combines the prediction results of these base learners and outputs [28]. In the construction of the RF model, two important parameters need to be set, namely the number of decision trees (ntree) and the number of variables randomly selected by the tree node (mtry). Generally, when the number of decision trees is greater than 500, the overall error rate of the model tends to be stable [29]. Therefore, the default value of ntree was 500, but it still needs to be further specified based on specific data. To ensure that the reliability of the prediction results did not affect the computational efficiency, the number of decision trees is set to 1000 in this study [30]. Generally, the number of variables randomly selected from decision tree nodes is one-third of the total number of independent variables [31]. However, the default value of mtry was not always able to be ascertained for the optimal model. It is necessary to tune the parameters [32]. Since 29 predictive variables were used in this study (1 ≤ mtry ≤ 29), a resampling technique, 10-fold cross-validation, was implemented to train and evaluate these 29 RF models (ntree = 1000, mtry = 1, 2, 3, …, 29). The model with the best generalization ability was selected as the optimal model. In this study, the randomForest package was used for RF modeling, and the caret package was used for 10-fold cross-validation in statistical software R4.0.3.

In this study, we did not test collinearity among these variables (Table A2). Although variable selection is an important part of machine learning, the tree models are resistant to uninformative predictors because the RF model selected in this study belongs to the decision tree model of ensemble learning. For example, if an independent variable is not used for any node splits during tree construction, the predictor model should not be functionally affected by that variable [32]. Also, RF models are insensitive to collinearity among independent variables, meaning the model’s results are resistant to unbalanced and missing data [31].

2.4.2. Model Validation

Stand age, topographic features, edaphic variables, and bioclimatic factors were used as independent variables, and volume gross growth was used as the dependent variable. The RF regression method established the volume gross growth model of each dominant tree species (group) of different origins. The regression model was evaluated using the following indicators [28]: coefficient of determination (R²), root mean square error (RMSE), relative root mean square error (rRMSE), and mean absolute error (MAE) of 10-fold cross-validation. The best-fit RF models for volume gross growth of each dominant tree species (group) in natural and planted forests showed various R² (0.370–0.922), RMSE (0.379–4.522), rRMSE (0.575–5.985), and MAE (23.22–55.81%) (

Table A4. Optimal model evaluation indicators of volume gross growth of dominant tree species (group) in natural and planted forests.

Origin	Dominant Species (Group)	Sample	Mtry	10-Fold Cross-Validation Results
				R2	MAE	RMSE	rRMSE
Natural forest	Pinus massoniana	2437	10	0.578	1.421	1.832	35.54%
	Pinus yunnanensis	1061	15	0.754	0.809	1.294	45.78%
	Pinus densata	221	26	0.449	1.330	1.734	46.21%
	Cunninghamia lanceolata	561	4	0.370	2.014	2.579	50.37%
	Cupressus	388	20	0.650	0.464	0.663	54.48%
	Quercus	984	29	0.661	0.506	0.726	30.81%
	Other hard broad leaf	68	12	0.473	1.415	1.796	48.21%
	Populus	40	27	0.603	1.172	1.422	40.07%
	Other soft broad leaf	105	7	0.464	1.243	1.603	55.81%
	Mixed coniferous	362	16	0.684	1.350	1.913	39.23%
	Mixed broad leaf forest	2766	22	0.922	0.379	0.575	23.22%
	Mixed coniferous and broad leaf forest	766	4	0.529	1.371	1.738	35.03%
Planted forest	Picea	67	23	0.528	1.671	2.101	51.71%
	Pinus armandi	169	8	0.481	2.247	2.932	44.82%
	Pinus massoniana	768	8	0.432	1.919	2.584	46.17%
	Pinus yunnanensis	159	11	0.519	1.343	1.737	38.93%
	Cunninghamia lanceolata	1203	8	0.380	1.010	1.246	26.73%
	Cryptomeria fortunei	125	21	0.417	4.522	5.985	55.13%
	Cupressus	402	29	0.634	0.964	1.288	39.12%
	Eucalyptus	98	1	0.468	2.511	2.520	44.62%
	Other hard broad leaf	148	10	0.614	2.100	2.741	46.85%
	Other soft broad leaf	85	7	0.504	1.399	1.774	41.41%
	Mixed coniferous	152	11	0.688	1.487	1.876	35.50%
	Mixed broad leaf forest	73	7	0.655	1.146	1.413	38.99%
	Mixed coniferous and broad leaf forest	105	3	0.405	1.465	1.854	32.60%

Notes: mtry, the optimized random forest parameter; R², squared correlation coefficient; MAE, mean absolute deviation; RMSE, root mean square error; rRMSE, relative root mean square error.

).

Similarly, stand age, dominant tree species, topographic features, edaphic variables, and bioclimatic factors were used as independent variables. The RF classification method was used to establish a loss judgment model and then determine whether there was a volume loss in each plot due to natural disturbances (e.g., withering, windfall, pests, and fire). The RF regression method established the arbor forest volume loss model. The classification model is evaluated as the following indicators [33]: True Positive, False Positive, True Negative, and False Negative (

Table A5. Confused matrix of random forest model.

Actual Type	Predictive Type
	Mortality	Non-Mortality
Mortality	True positive	False negative
Non-mortality	False positive	True negative

Notes: True Positive is the number of samples predicted by the model to be lost and the number of samples actually lost, False Positive is the number of samples predicted by the model to be lost but not actually lost, True Negative is the number of samples predicted by the model that no loss and no loss will actually occur, and False Negative is the number of samples that the model predicts will not be lost but actually are lost (Li et al., 2012).

). The accuracy index is calculated to test the model performance [34]. The prediction result of the RF classification model in natural and planted forests exhibited a relatively low misjudgment rate (0.36–4.63%) and then resulted in high model performance (R² = 0.712 and 0.763, respectively) in the RF model for volume gross loss for natural and planted forest (

Table A6. Mortality model results in dominant tree species (group) in natural and planted forests.

Origin	Actual Type	Predictive Type		Error Rate/%	Model Results
		Mortality	Non-Mortality		R2	MAE	RMSE	rRMSE
Natural forest	Mortality	6733	24	0.36%	0.712	0.544	0.974	95.20%
	Non-mortality	162	5235	3.00%
Planted forest	Mortality	3974	193	4.63%	0.763	0.870	1.407	84.77%
	Non-mortality	32	1387	2.26%

).

The net volume growth of arbor forests in the YREB has no significant difference between the predicted value by the RF model and the monitored value in the NFI (Figure 2). The net volume growth of arbor forests differed by 3.7% (between 7th and 8th NFI) and −8.9% (between 8th and 9th NFI) (p < 0.05), respectively. The difference between positive value and negative value was possibly related to edaphic variables. A previous study showed that soil nutrients (e.g., carbon, nitrogen, and phosphorus) varied over time across subtropical China [35]. But we assumed unchanged soil nutrients in this study. Overall, our results showed that the RF model established in this study could be used to predict the net volume growth of arbor forests in the YREB.

2.4.3. Prediction of Forest Volume Growth and Loss

This study only predicted future stand volume growth and volume loss of arbors in the existing forests and the new afforestation forests using the established RF model, based on the stand age, topographic features, edaphic variables, and bioclimatic factors. For arbor forest, the volume net growth in each permanent plot (Equation (1)), the future total volume in some year (Equation (2)), and the future total volume of each province in some year (Equation (3)) were estimated using following equations, respectively.

$$\Delta V\_NetGrowt{h}_{p,t}=\Delta V\_GrossGrowt{h}_{p,t}-\Delta V\_los{s}_{p,t}$$

(1)

$$V_{p,t}=V_{p,t-1}+\Delta V\_NetGrowt{h}_{p,t}$$

(2)

$$V_{i,t}={\displaystyle \sum }_p\left(V_{p,t}\times A_p\right)$$

(3)

where $\Delta V\_NetGrowt{h}_{p,t}$ is the volume net growth of the p-th plot in t year (m³ hm⁻² yr⁻¹), $\Delta V\_GrossGrowt{h}_{p,t}$ is the volume gross growth of the p-th plot in t year (m³ hm⁻² yr⁻¹), $\Delta V\_los{s}_{p,t}$ is the volume loss of the p-th plot in t year (m³ hm⁻² yr⁻¹), $V_{p,t}$ is the future total volume of the p-th plot in t year (m³ hm⁻²), $A_{p,t-1}$ is the total volume of the p-th plot in t-1 year (m³ hm⁻²), $V_{i,t}$ is the future total volume of arbor forest in t year in province i (m³), and $A_p$ is the grid area represented by p-th plot.

2.5. Estimation of Forest Biomass Carbon Storage and Annual Carbon Sink

The dominant tree species (groups) with their physiological similarities in the various plots of the NFI were combined (

Table A7. Biomass estimation model and carbon accounting factors for main arbor forests in YREB.

Dominant Species (Group)	Biomass Estimation Model B = (a·Vb)·λ					Carbon Fraction
	Sample Number	a	b	Correlation Coefficient R2	Correction Coefficient λ
Picea/Abies	25	5.413	0.633	0.966	1.012	0.493
Pinus massoniana	64	2.28	0.779	0.93	1.032	0.525
Cunninghamia lanceolata	199	4.012	0.631	0.924	1.018	0.506
Cupressus	26	6.711	0.569	0.758	1.04	0.51
Quercus	18	1.682	0.918	0.978	1.007	0.48
Other hard broad leaf	59	3.3	0.741	0.884	1.035	0.476
Populus	19	1.703	0.803	0.884	1.027	0.491
Eucalyptus	34	3.01	0.715	0.774	1.028	0.491
Other soft broad leaf	32	4.366	0.688	0.846	1.055	0.491
Mixed coniferous	11	6.699	0.538	0.808	1.012	0.502
Mixed broad leaf forest	20	1.526	0.908	0.898	1.028	0.479
Mixed coniferous and broad leaf forest	54	3.088	0.734	0.832	1.033	0.494

). Based on the relationship between the volume per unit area and the biomass of arbor forests in the China Forest Biomass Database [36], a fitting equation (Equation (4) and (5)) of the main dominant tree species (groups) in the YREB was established (Table A7). And the estimated biomass density of various plots of arbor forest was used to calculate the forest biomass C storage using the following equation (Equation (6)). The annual change in the biomass C storage of the arbor forest is regarded as the annual biomass C sink of the arbor forest (Equation (7)).

$$\ln \left(B_j\right)=\ln \left(a_j\right)+b_j\times \ln \left(V_j\right)$$

(4)

$$B_{j,p}=a_j\times V_{j,p}{}^{b_j}\times {\lambda }_j$$

(5)

$$C_j={{\displaystyle \sum }}_p\left(B_{j,p}\times C{F}_j\times A_p\right)$$

(6)

$$\Delta C_t={{\displaystyle \sum }}_j{C}_{j,t}-{{\displaystyle \sum }}_j{C}_{j,t-1}$$

(7)

where $B_j$ was the biomass density of tree species j (Mg hm⁻²), $V_j$ was the volume of tree species j per unit area (m³ hm⁻²), $a_j$ and $b_j$ were the coefficients of the biomass-volume correlation equation of tree species j; $B_{j,p}$ was the biomass density of plot p of tree species j (Mg hm⁻²), $V_{j,p}$ was the volume of plot p of tree species j per unit area (m³ hm⁻²), ${\lambda }_j$ was the biomass density correction coefficient of tree species j; $C_j$ was the biomass C storage of tree species j (Mg), $C{F}_j$ was the biomass C fraction of tree species j; $A_p$ was the grid area represented by sample plot p (hm²); $\Delta C_t$ was the biomass C sink of the arbor forest in year t (Mg C yr⁻¹).

For shrub and bamboo forests, we assume that their biomass density remains constant. The average biomass density of shrub forests in the YREB is set to 19.76 Mg hm^−2, and its C factor is 0.5 [37]. The average biomass density of Phyllostachys pubescens forest and other bamboo forests are set to 81.9 Mg hm⁻² and 53.1 Mg hm⁻², respectively, and their C factor is 0.47 [38].

3. Results

3.1. Forest Biomass C Storage in the YREB

The total forest biomass C storage in the YREB ranged from 3053.27 Tg C in 2015 to 6392.19 Tg C in 2060 (Figure 3). The arbor forests had the highest C storage between 2015 to 2060, accounting for 89.7–95.1% of total forest biomass C storage in the YREB, and the shrub forest C storage was comparable to the bamboo forest C storage (

Table A8. Variations in C storage and C density of existing forests and new afforestation forests from 2015 to 2060 in the YREB.

Type	Year	C Storage (Tg C)				C Density (Mg C hm−2)
		Arbor	Shrub	Bamboo	Total	Arbor	Shrub	Bamboo	Total
Existing forests	2015	2739.54	159.76	153.97	3053.27	38.53	10.67	34.96	33.75
	2020	3164.08	159.76	153.97	3477.81	44.49	10.67	34.96	38.44
	2030	3910.13	159.76	153.97	4223.86	54.54	10.67	34.96	46.69
	2040	4583.92	159.76	153.97	4897.65	63.69	10.67	34.96	54.13
	2050	5214.01	159.76	153.97	5527.74	72.35	10.67	34.96	61.10
	2060	5817.90	159.76	153.97	6131.63	80.68	10.67	34.96	67.77
New afforestation forests	2020	27.31	0.00	1.06	28.37	14.16	3.71	13.06	14.11
	2030	128.77	0.00	6.37	135.14	23.75	7.95	28.01	23.92
	2040	273.50	0.01	11.15	284.65	31.75	8.76	30.87	31.72
	2050	438.96	0.01	15.14	454.10	40.89	9.54	33.62	40.59
	2060	574.30	0.01	15.68	589.99	53.49	9.88	34.83	52.74
Total	2015	2739.54	159.76	153.97	3053.27	38.53	10.67	34.96	33.75
	2020	3191.39	159.76	155.03	3506.18	43.70	10.67	34.57	37.91
	2030	4038.90	159.76	160.34	4359.00	52.78	10.67	34.62	45.35
	2040	4857.42	159.77	165.12	5182.30	60.94	10.67	34.65	52.11
	2050	5652.97	159.77	169.11	5981.85	69.08	10.67	34.84	58.84
	2060	6392.19	159.77	169.66	6721.61	78.11	10.67	34.95	66.12

). From 2015 to 2060, the total forest biomass C storage was ranked in descending order: Yunnan, Sichuan, Jiangxi, Hunan, Guizhou, Hubei, Zhejiang, Chongqing, Anhui, Jiangsu, and Shanghai (Figure 3). The predicted forest biomass C storage in all provinces and municipalities linearly increased over time (Figure 3), but the proportions of total forest biomass C storage decreased in two provinces: from 35.6% in 2015 to 29.5% in 2060 for Yunnan province and from 26.4% in 2015 to 21.5% in 2060 for Sichuan province (Figure 4). Moreover, both current and future C storage in the existing forests made a greater contribution to the total forest biomass C storage in the YREB than those in the new afforestation forests (Table A8).

3.2. Forest Biomass C Density in the YREB

The total forest biomass C density in the YREB increased from 33.75 Mg C hm⁻² in 2015 to 66.12 Mg C hm⁻² in 2060 (Figure 5a), with an average of 49.01 ± 12.40 Mg C hm⁻² (Figure 5b). Overall, the predicted forest biomass C density in all provinces and municipalities significantly increased from 2015 to 2060 (Figure 5a), the top 2 regions with the highest increasing rate were the following: Zhejiang (from 28.27 to 76.53 Mg C hm⁻²) and Jiangxi (from 29.59 to 77.15 Mg C hm⁻²), while the top 2 regions with the lowest increasing rate were the following: Shanghai (from 25.90 to 46.87 Mg C hm⁻²) and Jiangsu (from 22.72 to 38.47 Mg C hm⁻²). The forest biomass C density varied significantly between different provinces and municipalities (p < 0.001), and the average forest biomass C density was ranked in descending order: Sichuan, Yunnan, Chongqing, Jiangxi, Zhejiang, Guizhou, Hubei, Anhui, Hunan, Shanghai, and Jiangsu (Figure 5b). Moreover, most provinces and municipalities had a higher average forest biomass C density than the YREB (Figure 5b).

3.3. Forest Biomass C Sink in the YREB

The forest biomass C sink in the YREB decreased from 90.58 Tg C yr⁻¹ in 2016–2020 to 73.98 Tg C yr⁻¹ in 2051–2060 (Figure 6a), with an average of 82.43 ± 6.17 Tg C yr⁻¹ (Figure 6b). Overall, the predicted forest biomass C sink in all provinces and municipalities linearly decreased from 2015 to 2060 (Figure 6a); the top 3 provinces with the lowest decreasing rate were the following: Shanghai (from 0.05 to 0.02 Tg C yr⁻¹), Chongqing (from 3.12 to 2.65 Tg C yr⁻¹), and Jiangsu (from 0.99 to 0.21 Tg C yr⁻¹), while the top 3 provinces with the highest decreasing rate were the following: Jiangxi (from 11.64 to 7.88 Tg C yr⁻¹), Hunan (from 10.65 to 7.32 Tg C yr⁻¹), and Hubei (from 7.69 to 6.06 Tg C yr⁻¹). The forest biomass C sink varied significantly between different provinces and municipalities (p < 0.001), and the average forest biomass C sink was ranked in descending order: Yunnan, Sichuan, Jiangxi, Hunan, Guizhou, Hubei, Zhejiang, Chongqing, Anhui, Jiangsu, and Shanghai (Figure 6b). Moreover, the C sink of new afforestation forests was predicted to reach the highest average of 16.95 Tg C yr⁻¹ in 2014–2050 (Figure A1).

4. Discussion

4.1. Forest Biomass C Sink Potential in the YREB

China’s key forestry ecological projects, dominated by afforestation, have increased land greening [17]. Therefore, afforestation is one of the most feasible and effective options for offsetting the greenhouse gas emissions and mitigating climate warming [9,12,14,16]. Since the implementation of the Grain for Green Program, the Yangtze River and Zhujiang River Shelter Forest Projects, and the Natural Forest Protection Project at the end of the 20th century [15,35], the forest area in the YREB increased continuously and was predicted to reach a value of 92.7 million hectares in 2060 (Table A1 and Table A3). In this study, the predicted forest biomass C storage in the YREB (including arbor, shrub, and bamboo in the existing and new afforestation forests) will increase by 3.67 Pg C from 2015 to 2060 (Figure 3 and Table A8), which is equivalent to 120% of the C storage of the existing forests in 2015. Meanwhile, the forest biomass C density will increase from 33.75 Mg C hm⁻² in 2015 to 66.12 Mg C hm⁻² in 2060 (Figure 5a and Table A8). These results indicated that the forests in the YREB have a large C sink potential in the future.

A previous study predicted that the total CO₂ emission in China will increase first from 2086 Tg C in 2016 to 2980 Tg C in 2027, and then decrease to 1370 Tg C in 2050 [39]. As a result, China is expected to emit 88,700 Tg C of CO₂ emission between 2016 and 2050. According to the predictions of this study, it is conservatively estimated that the forest in the YREB will absorb 4.1% of the CO₂ emission in China during the same period, which is equivalent to 68.3% of China’s forest biomass C sink potential [9,12]. Moreover, the proportion of forest biomass C sink in the YREB to CO₂ emission in China decreased first from 3.2% in 2016–2020 to 2.9% in 2021–2030, and then increased to 4.3% in 2041–2050 (Figure 7a), suggesting higher forest biomass C sequestration efficiency in the future. However, forest biomass C sink gradually decreased over time (Figure 6a), similar to previous studies on China’s forest vegetation [9,12,20]. A reason is related to the high percentage of mature and over-mature forests in the future [8]. Although this study predicted that forest biomass C sink in the new afforestation forests increased from 5.67 Tg C yr⁻¹ in 2016–2020 to 16.95 Tg C yr⁻¹ in 2041–2050 (Figure A1), there is a low total area of new afforestation forests (Table A3). To achieve a high C sequestration rate in the future, therefore, it is important to perform more efficient forest management (such as forest structure adjustment and forest quality improvement) in the YREB.

Interestingly, the predicted average forest biomass C storage between 2020 to 2050 in the YREB was closer to that in the same region reported by Qiu et al. [12] rather than by He et al. [9] and Gu et al. [22] (Figure 7b). One possible reason was the differences in the vegetation types. This study and Qiu et al. [12] collected data from the NFI in China. They predicted forest biomass C storage in similar vegetation types (arbor, shrub, and bamboo forests for this study versus arbor, shrub, bamboo forests, and nursery land for Qiu et al. [12]. While He et al. [9] collected data from the field survey of China’s forests that only included arbor forests such as deciduous broadleaf forest, deciduous needleleaf forest, evergreen broadleaf forest, evergreen needleleaf and needleleaf, and broadleaf mixed forest. And Gu et al. [22] predicted C storage in vegetation types, including forest, shrubland, grassland, and cropland. In these studies, the differences in vegetation types will result in different forest areas and affect the prediction of forest C storage. Moreover, compared to this study and Qiu et al. [12], the underestimated forest area reported by He et al. [9] and the overestimated vegetation area reported by Gu et al. [22] would lead to a high C density and a comparable C density in the YREB, respectively (Figure 7b). On the other hand, a previous study showed that global forest C balance was mainly regulated by soil nutrient availability [12], indicating the limitation of low soil nutrients across subtropical China to forest C sequestration [35]. Therefore, forest C storage based on climate-induced prediction [8,26] may differ from that based on multivariate (e.g., climatic factors, soil properties, stand attributes, and topographic features) prediction (e.g., this study and Qiu et al. [12]). Additionally, identifying the effects of biotic and abiotic factors on forest C storage across scales and ecosystems is beneficial to further understanding these differences [40,41].

4.2. Differences in the Forest Biomass C Sink Potential Amongst Regions

Between 2015 and 2060, Yunnan and Sichuan provinces would have the greatest forest biomass C sink with a total value of 0.99 and 0.66 Pg C, respectively (Figure 3), and the highest forest biomass C density with the average value of 64.81 and 66.83 Mg C hm⁻², respectively (Figure 5b). Benefited from the complex landform (Figure 1), the subalpine and alpine area of the western region in the YREB (such as Yunnan and Sichuan in this study) represents the second largest natural forest in China [8,22], which has diverse vegetation types and high forest productivity [42]. As a result, these forests are the most important C sink in YREB (Figure 6). However, the proportion of forest biomass C storage in these provinces to the total forest biomass C storage in the YREB during the same period decreased over time (Figure 4). These results indicated that the forest biomass C sequestration rate in these regions might weaken in the future, as illustrated by the decreasing forest biomass C sink (Figure 6a). This may be related to the weakening of annual C sink capacity caused by the gradual maturity of existing forests and new forests in the future and the decline in growth rate [43]. Therefore, it is necessary to focus on strengthening forest protection and management in the upper reaches of the YREB to maintain and enhance its forest biomass C storage function.

Consistent with the results of previous studies [9,12,22], the lower reaches, such as Shanghai municipality and Jiangsu province in the Yangtze River, had the lowest forest biomass C storage, C density, and C sink in the historical period and in the future (Figure 3, Figure 4, Figure 5 and Figure 6). Two mechanisms, a large proportion of young- and middle-age arbor forests and the low forest area, were the main reason for the lowest forest biomass C sink potential in these regions [8]. In addition, there is dense population, lack of forest resources, huge man-made destruction, and serious utilization interference, which will result in poor forest quality and low C sequestration. Surprisingly, forest biomass C density in these two regions has shown a steady increase from 2015 to 2060 (Figure 5a), indicating that the C sink function of forest biomass is continuously increasing after improved forest management.

In this study, the predicted forest biomass C sink potential of various regions (Figure 3, Figure 4, Figure 5 and Figure 6) differed from previous studies [9,12] due to differences in the estimation methods, the estimated forest types, and the setting of new afforestation scenarios. Moreover, the forest biomass C sink potential of various regions was unbalanced (Figure 3, Figure 4, Figure 5 and Figure 6). Overall, the highest forest biomass C density was in the upper reaches of the Yangtze River (e.g., Sichuan and Yunnan), followed by the middle reaches (e.g., Chongqing, Hubei, Hunan, Guizhou, and Jiangxi), and the lowest was in the lower reaches of the Yangtze (e.g., Anhui, Jiangsu, and Shanghai). These results were closely related to the forest coverage and the economic development level in these regions and showed an important function in increasing forest biomass C sinks and reducing emissions in the upper and middle reaches of the Yangtze River in the future.

4.3. Uncertainty and Limitation of This Study

In this study, although we predicted the increased forest biomass C sink potential in the YREB (Figure 7a), there was still some uncertainty. On the one hand, since future deforestation and logging will be greatly affected by the policy, it is difficult to predict and has high uncertainty. Therefore, the scenario set in this study is to assume that the forest will not be subjected to deforestation and logged in the future, and the forest will grow according to its natural growth process. Indeed, if the forests encounter the effects of anthropogenic and natural disturbances that lead to deforestation or death during the growth process, low-density young forests will replace high-density mature forests, which will cause deviations in the prediction results. On the other hand, one of the scenarios assumed that the soil fertility in all permanent plots remains constant in the future. Unfortunately, soil nutrients changed over time in subtropical China, resulting in uncertain predictions [13,35]. Moreover, the total area of new afforestation forests in each province and municipality and the allocation ratios of forest area in different periods, different origins, and different tree species may change with the adjustment of policies and management methods in the future, suggesting the changed C sink potential of new afforestation forests. Therefore, how to accurately estimate and predict forest biomass C sink potential in the YREB needs to be further studied. More importantly, taking into account the C sink potential of dead organic matter C pools (litter and deadwood), soil organic C pools, and harvested wood forest products C pools is necessary to assess regional carbon neutrality [22]. Furthermore, due to the lack of reliable C emission data in the YREB, we cannot make an accurate estimate of the future contribution of forest biomass C sink in offsetting CO₂ emissions in this region and only estimate the national contribution. In future research, it is necessary to strengthen the research on C emissions in different regions of the country.

5. Conclusions

Based on the NFI data from permanent plots, the RF model combined with multiple variables, such as regional climate, soil, stand, and topography, was used to predict the forest biomass C sink potential in China’s YREB from 2015 to 2060. The predicted forest biomass C storage and C density gradually increased over time, while the predicted forest biomass C sink gradually decreased over time in the future. These three C related indicators significantly differed across provinces and municipalities, the upper reaches (e.g., Sichuan and Yunnan), the middle reaches (e.g., Chongqing, Hubei, Hunan, Guizhou, and Jiangxi), and the lower reaches (e.g., Anhui, Jiangsu, and Shanghai) were listed in descending order. Overall, the C storage in the existing and new afforestation forests would increase by 3.67 Pg C from 2015 to 2060, equivalent to 120% of the C storage of the existing forests in 2015. It is conservatively estimated that the forest in the YREB will absorb 4.1% of the CO₂ emission in China between 2016 and 2050. These results indicated higher forest biomass C sequestration efficiency in the YREB in the future. Our results also suggested that improved forest management in the upper and middle reaches of the Yangtze River will help to enhance forest biomass C sink in the future.