1.1 study area
The Lijiang River Basin is located at 110°5'-110°44'E and 24°38'-25°56'N, with an area of about 5768.29km2, which is shown in Fig. 1. It originates from the south of the old boundary in the northeast of Maoer Mountain in Xing'an County, Guilin, and flows from north to south into Pingle Town and Gongcheng River in Pingle County, with a total length of 214 km [19]. The basin has a subtropical monsoon climate, with an annual average temperature of around 20°C and an annual average rainfall of 1949.5 mm. The four seasons are distinct, characterized by hot and rainy summers, and mild and rainy winters[20–22].
1.2 Methods
1.2.1 Land use simulation based on Markov-PLUS model
Compared with the traditional CA model, the PLUS model is based on the mining framework of land expansion analysis strategy (LEAS) while the CA model is based on multi-type random seeds (CARS). This model can mine the driving factors of land expansion and landscape change, and obtain higher simulation accuracy and more similar landscapes [23, 24] .
(1) land expansion analysis strategy
The random forest algorithm is used to mine the growth probability and driving factors to obtain the contribution of the development probability and driving factors of various types of land.
\(P_{{(i,k(x))}}^{d}= - P_{{(i,k(x))}}^{d}=\frac{{\sum\nolimits_{{n=1}}^{M} {I({h_n}(x)=d)} }}{M}\)
where \(P_{{(i,k(x))}}^{d}\)is the development probability of -type land in unit , and \(d=1\), if other land use types converted toland use types, \(d=0\)otherwise. is a vector composed of multiple driving factors; is the indicator function of the decision tree; is the total number of decision trees; \({h_n}(x)\) is the prediction type of the nth decision tree of vector .
(2) CA model based on multi-type random seeds
Under the constraints of development probability, combined with random seed generation and threshold decreasing mechanism, the PLUS model can automatically generate spatio-temporal dynamic simulation patches. Using the Monte Carlo method, when the land neighborhood effect is 0, the probability surface of each land \(OP_{{i,k}}^{{1,t}}\) as follows
\(OP_{{i,k}}^{{1,t}}=\left\{ \begin{gathered} P_{{i,k}}^{1} \times (r \times {\mu _k}) \times D_{k}^{t}\;,\;if\;\Omega _{{i,k}}^{t}=0\;and\;r<P_{{i,k}}^{1}, \hfill \\ P_{{i,k}}^{1} \times \;\Omega _{{i,k}}^{t} \times D_{k}^{t},\;all\;others. \hfill \\ \end{gathered} \right.\)
\(\Omega _{{i,k}}^{t}=\frac{{con(c_{i}^{{t - 1}}=k)}}{{n \times n - 1}} \times {w_k}\)
\(D_{k}^{t}=\left\{ \begin{gathered} D_{k}^{{t{\text{-1}}}},\;if\;\left| {G_{k}^{{t{\text{-1}}}}} \right| \leqslant \left| {G_{k}^{{t{\text{-2}}}}} \right|, \hfill \\ D_{k}^{{t{\text{-1}}}} \times \frac{{G_{k}^{{t{\text{-2}}}}}}{{G_{k}^{{t{\text{-1}}}}}},\;if\;0>G_{k}^{{t{\text{-2}}}}>G_{k}^{{t{\text{-1}}}},\; \hfill \\ D_{k}^{{t{\text{-1}}}} \times \frac{{G_{k}^{{t{\text{-2}}}}}}{{G_{k}^{{t{\text{-1}}}}}},\;if\;G_{k}^{{t{\text{-1}}}}>G_{k}^{{t{\text{-2}}}}>0. \hfill \\ \end{gathered} \right.\)
where is random values in the range of 0–1, \({\mu _k}\) is the threshold value of new land use patches generated by -type land use type, \(D_{k}^{t}\) is an adaptive drive coefficient and it depends on the gap between the current number of land use type and the target demand, \(\Omega _{{i,k}}^{t}\) is the neighborhood effect of the unit which means that \(\Omega _{{i,k}}^{t}\) is the component of the neighborhood under the proportion of land cover.is the type of land use.
1.2.2 Carbon storage calculation based on InVEST model
Based on the landscape type data and carbon density table (see Table 1) of the Lijiang River Basin, the carbon storage and sequestration module in the InVEST model [25,26] was used to evaluate the change of carbon storage in the Lijiang River Basin, the function are as follows
Ci_tot=Ci_above+Ci_below+Ci_soil+Ci_dead (1)
where is the land use type, \({C_{i\_tot}}\) is the total carbon storage in the ecosystem; \({C_{i\_above}}\) is the aboveground carbon storage of vegetation(such as bark, trunk, leaves, etc.); \({C_{i\_below}}\) is the underground carbon storage of vegetation (vegetation roots); \({C_{i\_soil}}\) is soil carbon storage (including mineral soil and organic soil); Ci _ dead is the carbon stock of dead organic matter (including dead trees and litter). Multiplying the biomass by carbon content of plant and combining with the vegetation biomass distribution map of China (spatial resolution of 500m), the carbon storage of plant is obtained. The carbon content of plant and soil was determined by concentrated sulfuric acid-potassium dichromate hydration heating method. Soil carbon storage (\({C_{soil}}\)) is calculated according to Teng et al[27]. Other parts of carbon storage refer to the relevant literature [28–30].
Table 1
Carbon density of land use landscape types in Lijiang River(t/hm2)
types
|
\({C_{i\_above}}\)
|
\({C_{i\_below}}\)
|
\({C_{i\_soil}}\)
|
\({C_{i\_dead}}\)
|
cropland
|
13.50
|
2.70
|
96.59
|
1.00
|
forest
|
73.59
|
25.20
|
207.33
|
3.50
|
shrub
|
18.96
|
5.69
|
9.4
|
2.47
|
grassland
|
5.01
|
13.53
|
117.06
|
1.00
|
waters
|
0.21
|
0.00
|
0.00
|
0.00
|
Barren land
|
19.52
|
3.9
|
0.86
|
0.00
|
impervious land
|
1.20
|
0.93
|
12.48
|
0.00
|
1.2.3 Simulation of future scenario
In order to ensure the sustainability of carbon storage resources in the Lijiang River Basin. Inertial development model, ecological priority model and urban development model are considered.
(1) Inertial scenario. Given the assumption that the evolution of landscape patterns in the Basin is unaffected by the new policy, Besides, the parameters for neighborhood factors and transfer cost matrices of each land use type remain unchanged.
(2) Ecological protection scenario. Emphasis has been placed on strengthening the management of ecological protection areas such as forests, grasslands, and water bodies. Measures have been implemented to curb the transformation of park green spaces and land within ecological protection areas into impervious land. There is a strict control in place to counter the disorderly expansion trend of impervious land, with the aim of mitigating the decline in forest and arable land. Therefore, in this scenario, decrease the conversion rate of forest, shrub, and grassland by approximately 50% to impervious land and increase the conversion rate of arable land to forest by 30% based on the inertial scenario and the land-use transfer matrix for the period 2021–2031.
(3) Urban development scenario. " The New Urbanization Plan for Guangxi (2021–2035)" indicates that the level of urbanization development in Guangxi is currently below the national average. There is a need for further improvement in urban and rural infrastructure to enhance the overall urbanization process. In this scenario, raise the conversion rate of cropland, forest, shrub, and grassland to impervious land by 20%. Besides, reduce the conversion rate of impervious land to arable land, forest, shrub, and grassland by 20%.
1.2.4 Data source
The land cover dataset for the Lijiang River Basin is derived from the annual land cover dataset published by Huang[31]. The factors influencing changes in land types mainly include natural and socio-economic factors, with specific indicators and sources outlined in Table 2.
Table 2
Data type
|
Name
|
Years
|
Source
|
Resolution
|
Topography and Landforms
|
DEM
|
2006
|
ALOS
|
12.5m
|
Slope
|
From DEM
|
Aspect
|
Curvature
|
Climate
|
Surface runoff
|
2010
|
ERA5
|
11132m
|
Temperature
|
Precipitation
|
ET
|
MOD16A2
|
500m
|
Distance
|
Distance to river
|
2015
|
National Geomatics Center of China
|
30m
|
Distance to artificial
|
Distance to town
|
2010
|
Distance to road
|
Distance to rural residential area
|
Vegetation
|
LAI
|
2010
|
MODIS
Landsat8 OLI/TIRS
|
500m
|
GPP
|
VCF
|
30m
|
Population Income
|
population density
|
2010
|
WoridPop
|
100m
|
Rural per capita net income
|
Social and Economic Statistical Yearbook
|
1000m
|
Fiscal Revenue
|
Agricultural development
|
Grain yield per unit area
|
\
|
Meat production per unit area
|
Tourism
|
Tourist point density
|
2020
|
Qunar combines Python and ArcGIS for kernel density calculation.
|
\
|
Hotel density
|
Else
|
LISA-I
|
2010
|
Land use type data calculation
|
\
|
Power plant density
|
2018
|
WRI kernel density calculation
|
\
|
Rocky desertification
|
2010
|
Combined with slope climate vegetation and other factors to calculate
|
30m
|