Spatial-temporal distribution pattern and driving factors of agricultural carbon sinks in Beijing-Tianjin-Hebei region from the perspective of carbon neutrality

Abstract

Developing policies to reduce carbon emissions in agriculture is crucial for achieving the ‘dual carbon’ goal. Therefore, scientifically analysing the temporal and spatial distribution characteristics of agricultural carbon sinks and their driving factors holds paramount importance for the coordinated and integrated development of regional agriculture and the realization of sustainable development. Based on the perspective of carbon cycle in agricultural production, the measurement system of agricultural net carbon sink was established from the perspective of carbon sink/carbon source, this study conducted an analysis of the temporal and spatial variation characteristics as well as driving mechanisms of agricultural net carbon sinks. The findings are as follows: (1) The agricultural net carbon sink exhibited an increasing trend from 2009 to 2020, with favourable intensity and level. (2) Agricultural land use and livestock and poultry production constituted the primary sources of agricultural carbon emissions. Notably, agricultural carbon emissions demonstrated a decreasing trend during the study period. (3) The net agricultural carbon sink displayed local spatial aggregation, indicating significant regional differences. (4) Factors such as rural economic level, urbanization level and agricultural employment significantly promoted the net carbon sink. In contrast, the rural industrial structure, technical level and crop seeding had inhibitory effects. Therefore, it is imperative to promote the reduction of agricultural carbon emissions in the Beijing-Tianjin-Hebei region. This entails accelerating the construction of new agriculture and rural areas, facilitating industrial upgrading, promoting the development of low-carbon-sink regions into high-carbon-sink regions and actively fostering the coordinated and integrated development of regional agriculture.

Introduction

In recent years, the rapid surge in greenhouse gas (GHG) emissions and the intensification of the global greenhouse effect have rendered climate change a pivotal threat to the survival and development of humanity. According to the Sixth Assessment Report of the Intergovernmental Panel on Climate Change, global surface temperatures in 2011–2020 were 1.1°C higher than in 1850–1900 (IPCC, 2022), with human activities very likely serving as the primary driver of global warming. Since the inception of the Kyoto Protocol, carbon emissions have garnered significant attention from the international community and academia (Cui et al., 2021). In this context, measures such as energy conservation, emission reduction and carbon sink compensation have emerged as pivotal strategies for achieving carbon neutrality (Zhang and Deng, 2022). In 2020, China established the goals of ‘peaking carbon’ by 2030 and achieving ‘carbon neutrality’ by 2060, aiming to harmonize economic development with environmental protection and facilitate the seamless realization of the ‘dual carbon’ goal (Shi et al., 2022). Currently, China finds itself in a critical phase transitioning from traditional and primarily modernized agriculture to comprehensively modernized agriculture (Tian and Chen, 2021). While the utilization of fertilizers, pesticides and agricultural machinery has boosted crop yields, the rapid increase in agricultural production efficiency has concurrently led to a significant rise in carbon emissions (Wang et al., 2020). To address the environmental and climate challenges arising from the swift development of agriculture and the associated increase in carbon emissions, the implementation of low-carbon emission reduction measures in the agricultural production process (Yu and Mao, 2022), and the enhancement of its carbon sink effect (Zhang and He 2022), have become primary strategies to combat global climate change and environmental crises (Liu et al., 2022). Therefore, exercising control over the carbon emissions of agricultural production factors, particularly harnessing the carbon sink effect of agricultural ecosystems (Li et al., 2021), will facilitate the promotion of low-carbon, green and healthy development in agricultural production, playing a pivotal role in advancing the early realization of the national ‘double carbon’ goal (Luo, 2016).

Agricultural carbon emission refers to the direct or indirect GHG emissions caused by the use of agricultural production factors and waste disposal in the process of agricultural production (He et al., 2018). As the world's second-largest source of carbon emissions (Jia et al., 2011), according to the IPCC report, agriculture and food systems were responsible for approximately 21–37% of total GHG emissions during the period 2007–2016 (Valérie et al., 2019). Agriculture production process has the dual attributes of being a carbon sink and carbon source (Chen et al., 2015a), and improving the agricultural carbon sink capacity and adjusting the carbon sink and carbon source structure have become major measures by which to mitigate agricultural carbon emissions (Cui et al., 2022b). At present, the main research directions in the field of agricultural carbon sinks are the spatiotemporal evolution characteristics of agricultural carbon sinks/sources (Xie and Liu, 2022) and the analysis of the driving factors of agricultural carbon sinks (Fang et al., 2007; Kuang et al., 2010). The net primary productivity model (Tian and Zhang, 2013), material balance method (Han et al., 2018), emission factor method (Yan et al., 2018) and direct measurement method have been used to estimate agricultural carbon sinks/sources. Among them, direct measurement method and material balance method are based on specific facilities and experimental processes to calculate carbon emissions, which have high measurement accuracy and can reflect the real carbon emissions, but the process is complex and not suitable for larger scale research (Liu et al., 2017; Qiu et al., 2022). When calculating agricultural carbon sink/carbon source at regional scale, emission factor method based on carbon source classification and net primary productivity model based on crop yield are mainly used at present. Several scholars examined the changes and spatiotemporal characteristics of agricultural carbon emissions in China over the 20 years following 1997 at the provincial scale (Huang et al., 2019). Shan et al. (2022) analysed the carbon emission efficiency and its influencing factors at the municipal scale in Hubei Province. Temporal and spatial evolution is an important means to study agricultural carbon sink/carbon source, and to deal with the possible ‘carbon crisis’ by understanding its changing characteristics and trends. Current studies on spatial and temporal change mainly focus on the analysis of the evolution trend of single carbon convergence at the national and provincial scales (Jiang et al., 2010; Cao et al., 2018, 2022). It mainly shows the characteristics of regional carbon sink/carbon source changes over timescale. Li et al. (2022b) studied the change characteristics of agricultural production efficiency and carbon sink in China from 2000 to 2019 by using data envelopment method. The analysis of driving factors can decompose the driving mechanism of temporal and spatial changes of agricultural carbon sink, which is of great significance for formulating coping strategies and means. At present, the main methods used are the Spatial Durbin Model (SDM), the Kaya identity model, Data Envelopment Analysis (DEA)-Malmquist and other research methods (Liu and Zhang, 2020; Guo et al., 2021; Zhu and Huo, 2022). For example, Huang and Zhu (2022) and Tian and Zhang (2020) used Geographically Weighted Regression Model (GTWR) and the log-average D-exponential decomposition method to analyse the influencing factors of the agricultural carbon effect. Liu et al. (2022) and Liu and Gao (2022) used the minimum distance to weak efficient frontier model to calculate the carbon emissions in the agricultural areas of the Yangtze River Economic Corridor, and conducted an empirical analysis on the influencing factors of agricultural carbon emissions by constructing the Tobit model. Deng et al. (2023) applied stochastic forest algorithm to quantitative analysis of factors affecting agricultural carbon emission efficiency, and achieved good results. In summary, the research method based on the decomposition of influencing factors is the most widely used method at present.

Several scholars have focused largely on the study of agricultural carbon sinks, while there have been few studies on the agricultural net carbon sink capacity (Yu et al., 2022; Li et al., 2022a). Moreover, the research scale is mostly concentrated at the provincial level, and research at the urban-regional scale is lacking (Zheng et al., 2022). Net carbon sinks are an important premise for reflecting the carbon effect of agricultural production. Furthermore, the measurement of the carbon emissions of agricultural production activities is mainly focused on the carbon emissions of crop production; there is no consideration of the carbon emissions of other links of agricultural production activities, such as animal husbandry, straw burning and other agricultural production processes. On the issue of regional agricultural carbon emissions, it is not advisable to ignore the carbon absorption effect. Therefore, on the basis of the research on agricultural carbon emissions, it is of scientific significance to supplement the research on the capacity and level of carbon sinks and net carbon sinks at the municipal level, and understand their temporal and spatial distribution effects and development characteristics for exploring the regional agricultural development level and mitigation effects.

As one of the most open and populous urban agglomerations in China, the Beijing-Tianjin-Hebei (BTH) region is an important engine driving China's economic development. However, it also faces serious problems such as environmental pollution and resource shortage, and excessive carbon emissions, which restrict the speed of development. Moreover, as an indispensable grain producing area in the Bohai Rim Economic Belt, agricultural industry activities are intensive, and agricultural production is a critical source of carbon emissions. Therefore, the improvement of regional carbon sink capacity, dynamic monitoring and assessment of its carbon sink capacity and level can provide support for the integrated and coordinated development of the BTH region and the implementation of green and low-carbon industrial economy. Based on this, in the present work, the statistical data of agricultural production in the BTH region were combined, and the issue of agricultural carbon emissions was fully considered. Based on the Theil index and spatial autocorrelation analysis, the spatial and temporal evolution patterns of the net carbon sink capacity of agricultural production in the BTH region were analysed, and the Stochastic Impacts by Regression on Population, Affluence, and Technology (STIRPAT) model was extended to decompose the driving factors of the net carbon sink. Some development rules of the net carbon sink in the BTH region are discussed and proposed to promote the green, low-carbon, and sustainable development of regional agricultural production.

Materials and methods

Profile of the study area

The BTH region is situated in the Bohai Sea Economic Belt in northern China, covering a geographic range from 113°27′E to 119°50′E and 36°05′N to 42°40′N (Fig. 1). It encompasses two municipalities directly under the central government's control, namely Beijing and Tianjin, along with 11 prefecture-level cities in Hebei Province. The total area is approximately 208 000 km², constituting about 2% of China's total land area. The terrain is higher in the northwest and lower in the southeast, with the highest elevation reaching 2691 m. The landform is characterized by complexity and diversity, experiencing four distinct seasons and favourable photothermal conditions. The region features a typical warm temperate sub-humid continental monsoon climate, with an annual average precipitation of around 500 mm. In 2020, the Gross Domestic Product (GDP) of the BTH region reached 8.6 trillion RMB. Within this, the total agricultural output amounted to 700 billion RMB, reflecting a year-on-year increase of 10.7%. The regional grain output contributed to 5.8% of the national total.

Figure 1. Map of geographical location of the study area.

Source of data

The study covered a 12-year period from 2009 to 2020. Statistical data on crop yield, gross agricultural product and planted area in each region were obtained from the Hebei Agricultural and Rural Statistical Yearbook (Hebei provincial Bureau of Statistics, 2021), the Beijing Statistical Yearbook (Beijing Municipal Bureau of Statistics, 2021) and the Tianjin Statistical Yearbook (Tianjin Bureau of Statistics, 2021). The charts in this paper were analysed using ArcGIS 10.2, Origin, Excel, SPSS and other software.

Measuring methods of agricultural net carbon sink

Measurement of agricultural carbon emissions

Due to the interaction and circularity of agricultural production activities, the carbon cycle process becomes relatively complex. When cities are chosen as the research scale, carbon emission sources become numerous and challenging to uniformly quantify. To simplify the calculation process, the Ran et al. (2017) estimation method is employed, where agricultural carbon emissions are calculated as the sum of each carbon emission source multiplied by its corresponding carbon emission coefficient. As shown in Eqn (1):

(1)$$C = \sum {C_i} = \sum {T_i} \times \theta _i$$

where C is the total agricultural carbon emissions, C _i is the agricultural carbon emissions of carbon source class i in year t, T _i is the amount of carbon source class i in year t and θ _i is the carbon emission coefficient of the corresponding carbon source.

According to the characteristics of agricultural development in the BTH region and the timeliness and availability of agricultural data, agricultural carbon emissions were divided into the following four aspects: (1) carbon emissions caused by agricultural land production activities (Li et al., 2011); (2) carbon emissions from crop life activities (Duan et al., 2011); (3) carbon emissions generated by livestock and poultry breeding (Min and Hu 2012); (4) carbon emissions from the open burning of crop straw (Cao et al., 2005).

1. 1) Carbon emissions from agricultural land-use activities

Carbon emissions from agricultural land use refer to GHG emissions directly or indirectly caused by various production factors of agricultural activities. These emissions primarily include: (1) carbon emissions generated directly or indirectly by fertilizers and pesticides during soil decomposition; (2) indirect carbon emissions resulting from the consumption of electricity during agricultural irrigation activities; (3) carbon emissions produced during the use of agricultural materials (plastic film, diesel, etc.). The carbon emission coefficient for each production factor is provided in Table 1.

1. 2) Carbon emissions from crop life activities

GHG emissions from agricultural crops mainly include CH₄ emissions from paddy fields and N₂O emissions from wheat, corn and other crops. Due to the significant impact of climate and temperature differences on the growth habit of rice, there are evident regional variations in rice varieties and planting seasons. According to the research by Min and Hu (2012), the rice varieties of the BTH region are mainly mid-season varieties (single-season late rice, winter paddy rice and stubble rice). Combined with the greenhouse emission coefficients of major crops measured by domestic scholars via experimentation, to achieve uniform quantification, when CH₄ and N₂O are uniformly replaced by standard C, according to the Sixth Assessment Report of the IPCC, the greenhouse effects caused by 1 t CH₄ and N₂O are equivalent to 6.8182 and 81.2727 t C, respectively (Table 2).

1. 3) Carbon emissions from livestock and poultry farming

During livestock breeding and the life activities of livestock and poultry, intestinal fermentation will produce CH₄, while faecal discharge will produce CH₄ and N₂O emissions. The carbon emission calculation formula is shown in Eqns (2) and (3):

(2)$$CH_{4live} = \sum\limits_{i = 1}^n {T_i \times \phi _i} $$

where CH _4live refers to CH₄ emission from livestock breeding, T _i refers to the annual average feeding quantity of livestock and poultry of i species and ϕ _i refers to the CH₄ emission coefficient of i type of livestock.

(3)$$N_2O_{live} = \sum\limits_{i = 1}^n {T_i \times \varphi _i} $$

where N ₂O _live is the N₂O emission of livestock and poultry breeding, T _i is the annual average feeding quantity of the i type of livestock and poultry and φ_i is the N₂O emission coefficient of i type livestock and poultry; according to the main livestock breeding types in the BTH region, which mainly include cattle, sheep, pigs, rabbits and poultry (Table 3).

Table 1. Carbon emission coefficient of production factors

C, carbon.

Table 2. Greenhouse gas emission coefficients of crops

N₂O, nitrous oxide; CH₄, methane.

Table 3. Carbon emission coefficients^a of livestock and poultry breeding

^a The coefficients refer to the results of Min and Hu (2012) and Shang et al. (2015).

Due to the different feeding cycles of livestock and poultry, the average annual feeding quantity of livestock and poultry should be adjusted. Based on the calculation method proposed by Hu and Wang (2010), the annual average feeding quantity of livestock and poultry was adjusted based on the feedlot rate. For pigs, sheep, rabbits and poultry, for which the feedlot rate is ⩾1, the annual average feeding quantity was adjusted according to the feedlot amount, refer to Eqn (4) for details:

(4)$$T_i = Days_-alive_i \times \displaystyle{{N_i} \over {365}}$$

where T _i is the annual average breeding quantity of livestock or poultry type i, Days ₋alive _i is the average life cycle of livestock or poultry type i and N _i is the average annual output of livestock or poultry type i. The average life cycles of pigs, sheep, rabbits and poultry were considered to be 200, 210, 105 and 55 days, respectively.

For livestock whose exit rate is less than 1 (the exit rate of cattle was considered to be less than 1 in this study), the average annual breeding quantity was adjusted according to the stock quantity at the end of the year, as shown in Eqn (5):

(5)$$T_i = ( {L_{it} + L_{i( t-1) }} ) /2$$

where T _i is the annual average breeding quantity of livestock type i, L _it and L _i(t−1), respectively, represent the stock quantities of livestock type i and poultry at the end of the year and the stock quantity at the end of the year t − 1.

1. 4) Carbon emissions from straw burning

The process of straw incineration produces large amounts of GHG such as CO and CO₂. Due to the different proportions and yields of straw incineration, the estimation formula for the carbon emissions of straw incineration is established, as shown in Eqn (6):

(6)$$G = \sum M _i \times R \times W_i = \sum {V_i \times S_i \times R \times W} _i$$

where G is the carbon emissions of straw incineration, M _i is the straw yield of crop type i, R is the proportional coefficient of straw incineration, W _i is the carbon emission coefficient of the straw burning of crop type i and V _i is the yield of crop type i. Moreover, S _i is the ratio of the coefficient of the grain to the grass of crop type i, i.e. the ratio of the crop straw yield to the crop economic yield. According to the research of Cao et al. (2005), the proportion of open straw incineration in Hebei is 30%, while those of Beijing and Tianjin are 0. The carbon emission coefficient of straw burning is shown in Table 4.

Table 4. Carbon emission coefficients^a of crop straw burning

^a The carbon emission coefficients are based on the research by Liu et al. (2011).

C, carbon.

Measurement of the agricultural carbon sink

The agricultural carbon sink specifically refers to the carbon uptake during the entire life cycle of crop growth, namely the net primary production formed by crops through photosynthesis, i.e. biological yield. In this study, the method for estimating the carbon sink, based on the crop productivity model established by Tian and Zhang (2013), is shown in Eqn (7):

(7)$$C_t = \sum {C_i} = \sum {c_i} \times D_i = \sum {c_i \times Y_i \times ( {1-r_i} ) /H_i} $$

where C _t is the carbon uptake of crops, i.e. the carbon sink, C _i is the carbon sink of the crop type i and c _i, r _i and H _i, respectively, represent the carbon content coefficient, water-cut ratio and economic coefficient of crop type i. Moreover, D _i and Y _i are the total biomass and economic yield of crop type i, as shown in Table 5.

Table 5. Carbon content coefficient, economic coefficient and water-cut ratio of crops

‘Carbon content coefficient’ represents the ratio of carbon mass to crop yield. ‘Coefficient of economy’ signifies the ratio of agricultural economic yield to biomass yield. ‘Water-cut ratio’ denotes the ratio of water content to crop yield. All the above variables are ratios without dimensions. The coefficient mainly refers to the studies of Cao et al. (2018) and Han et al. (2012).

Measurement of the agricultural net carbon sink

The agricultural net carbon sink refers to the difference between the carbon sink and carbon emissions in agricultural production activities. To further analyse the temporal and spatial variation characteristics of the net carbon sink and the level of carbon sink intensity in the BTH region, the net agricultural carbon sink and cultivated land area were used to reflect the intensity of the agricultural net carbon sink. Additionally, the ratio of the carbon sink to total carbon emissions was used to reflect the carbon sink level (Lyu, 2019). The calculation formula is shown in Eqns (8)–(10):

(8)$$C_{net} = C_t-C_a$$

(9)$$C_s = C_{net}/S_{land}$$

(10)$$C_l = C_t/C_a$$

where C _net is the net agricultural carbon sink; C _t is the total amount of carbon absorbed by crops, also known as carbon sink; C _a is total agricultural carbon emission; C _s is the agricultural net carbon sink intensity; S _land is cultivated land area; C _l is the level of agricultural carbon sink.

Theil index

The Theil index is a special form of a generalized entropy index system used to measure income inequality between individuals or regions (Yuan and Liu, 2018); it demonstrates good decomposability and the ability to independently measure the contributions of intra-group and inter-group differences to total differences based on the ‘entropy’ theory in information theory. Its value is generally within the interval [0,1]; the smaller the value, the smaller the regional difference. Therefore, based on the original model improved by Cui et al. (2022a), the new Theil index and its decomposition model were obtained as shown in Eqns (11)–(18):

(11)$$T = \sum\limits_i {\left({\displaystyle{{C_i} \over C}} \right)} \ln \left({\displaystyle{{C_i/C} \over {X_i/X}}} \right)$$

(12)$$T_{wi} = \sum\limits_i {\left({\displaystyle{{C_{\,ji}} \over {{\rm C}_j}}} \right)} \ln \left({\displaystyle{{C_{\,ji}/C_j} \over {X_{\,ji}/X_j}}} \right)$$

(13)$$T_w = \sum\limits_j {\left({\displaystyle{{C_j} \over C}} \right)} T_{wi} = \sum\limits_j {\sum\limits_i {\left({\displaystyle{{C_j} \over C}} \right)} } \left({\displaystyle{{C_{\,ji}} \over {C_j}}} \right)\ln \left({\displaystyle{{C_{\,ji}/C_j} \over {X_{\,ji}/X_j}}} \right)$$

(14)$$T_b = \sum\limits_j {\left({\displaystyle{{C_j} \over C}} \right)} \ln \left({\displaystyle{{C_{\rm j}/C} \over {X_j/X}}} \right)$$

(15)$$T = T_w + T_b = \sum\limits_j {\left({\displaystyle{{C_j} \over C}} \right)} T_{wi} + T_b$$

(16)$$R_a = \displaystyle{{T_w} \over T}$$

(17)$$R_b = \displaystyle{{T_b} \over T}$$

(18)$$R_j = \displaystyle{{X_j} \over X}\displaystyle{{T_{wi}} \over T}$$

In these equations, T, T_wi, T_w and T_b represent the total unit of carbon sink Theil index (%), the carbon sink Theil index (%), the carbon sink Theil index (%) within the region and the carbon sink Theil index (%) between the regions, which can reflect the differences of the BTH region as a whole, each sub-region, within the region and between regions, respectively.

Moreover, C_i, C, C_ji and C_j are carbon sinks (t) by urban area, total carbon sinks (t) by BTH, carbon sinks by urban area (t) by sub-region and total carbon sinks by sub-region (t). X_i, X, X_ji and X_j represent the total agricultural output value (100 million RMB) or rural population (10 000 people) of each urban area, BTH region, city and sub-region, respectively. Eqn (15) indicates that the total difference is made up of intra-regional and inter-regional differences (%). R_a, R_b and R_j represent intra-regional, inter-regional and sub-regional contribution rates (%), respectively. When X represents population size, it represents the per capita carbon sink Theil index, expressed by T(P); when X represents GDP, it represents the carbon sink intensity Theil index, expressed by T(G).

Spatial autocorrelation analysis

Spatial autocorrelation analysis is used to measure whether phenomena exhibit agglomeration, dispersion or random distributions in space. In this study, the Global Moran's I Index was employed to test for correlations between net carbon sinks in the BTH region. The calculation formula is shown in Eqn (19):

(19)$$I = \displaystyle{{\sum\nolimits_{i = 1}^n {\sum\nolimits_{\,j = 1}^n {w_{ij}} } ( {x_i-\bar{x}} ) ( {x_j-\bar{x}} ) } \over {S^2\sum\nolimits_{i = 1}^n {\sum\nolimits_{\,j = 1}^n {w_{ij}} } }}$$

where I is the Global Moran's I Index, the range of which is [−1,1]. The closer the value is to 1, the higher the spatial positive correlation, and the closer the inter-regional connection. The closer the value is to −1, the higher the spatial negative correlation, and the further the inter-regional connection. A value of 0 indicates no spatial correlation. Moreover, n represents the number of regions, x_i and x_j represent the net carbon sink values of regions i and j, respectively; $\bar{x}$ represents the sample mean, S ² represents the attribute variance and w_ij represents the spatial weight matrix of regions i and j, constructed based on the inverse distance weight standard.

Local spatial autocorrelation can be used to explore whether similar or dissimilar index values are clustered together in a local region. Therefore, the Local Moran's I Index was used to analyse the agglomeration situation of the agricultural net carbon sink in the BTH region, and the calculation formula is shown in Eqn (20):

(20)$$I_i = \displaystyle{{( {x_i-\bar{x}} ) } \over {S^2}}\sum\limits_{\,j = 1}^n {w_{ij}} ( {x_j-\bar{x}} ) $$

where I_i is the Local Moran's I Index and the other letters have the same meaning as in Eqn (19) above. Positive values indicate the existence of high-high (H-H) or low-low (L-L) spatial agglomeration in this region, while negative values indicate the existence of high-low (H-L) or low-high (L-H) spatial agglomeration in this region. The greater the absolute value of the Local Moran's I Index, the higher the degree of spatial agglomeration.

STIRPAT model

The STIRPAT model, developed on the basis of the Impact, Population, Affluence, and Technology (IPAT) model, overcomes the quantity limit of the IPAT model and avoids the influence of the same scale change on the IPAT model (Huang et al., 2021). It is widely used to study the impact of population, economy and technology on the environment, reflecting the impact of economic and social development (Liu et al., 2023). The model has been studied by many scholars on carbon emissions and related issues (Jiang et al., 2023). The factors influencing the net carbon sink in the BTH region were examined using the extended STIRPAT model, with the standard form referred to Eqn (21):

(21)$$I = aP^bA^cT^de$$

where a is the model coefficient, and I, P, A and T are the environmental impact, population size, economic affluence and technological level, respectively. Moreover, ^b, ^c and ^d are the elastic coefficients of population, wealth and technology, and e is the error term. By taking the natural logarithms of both sides, the model is shown in Eqn (22):

(22)$$\ln ( I) = \ln ( a) + b \times \ln ( P) + c \times \ln ( A) + d \times \ln ( T) + \ln ( e) $$

where b, c and d are the regression coefficients of the corresponding explanatory variables, which reflect the influence of the changes of each driving factor on the dependent variable. Each 1% change in P, A and T will result in a b%, c% and d% change in I, respectively, while the other coefficients remain constant.

To further study the influences of factors such as rural economic development, agricultural technology level, agricultural industry structure, urbanization degree and rural population on the carbon sink capacity of the BTH region, based on existing research findings and the availability of data, we have selected eight indicators to expand the model. The model construction is shown in Eqn (23):

(23)$$CE = a \times P^b \times A^c \times T^d \times V^e \times U^f \times C^g \times R^h \times Y^i \times {\rm e}$$

By taking the logarithms of both sides of the equation, the formula can be obtained as shown in Eqn (24):

(24)$$\eqalign{\ln CE = & \ln a + b\ln P + c\ln A + d\ln T + e\ln V + f\ln U \cr & \quad + g\ln C + h\ln R + i\ln Y + \ln e} $$

where CE is the net agricultural carbon sink in the BTH region (10⁴ t), P is the total agricultural output value (AOV) (100 million RMB) and A is the proportion of the agricultural output value among the total output value of agriculture, forestry, husbandry and fishery (AP) (%). Moreover, T is the urbanization rate, which is expressed as the proportion of the urban population among the total population of the BTH region (UR) (%), and V is the technical factor, which represents the development level of agricultural technology and is represented by the total power of agricultural machinery (AMP) (10⁴ kW). Additionally, U is the sown area of crops (ACR) (hm²), C is the number of agricultural employees (AE) (10 000 people), R indicates the development level of the agricultural economy, expressed by the per capita agricultural output value of rural areas (ARV) (RMB/person), and Y represents the affected area of farmland (AFA) (hm²).

The multi-collinearity between factors is the primary shortcoming of the STIRPAT model, addressed through partial least-squares regression. Utilizing SPSS software, dimensionality reduction analysis extracted the principal components of eight influencing factors: total agricultural output value, proportion of the agricultural output value, urbanization rate, total agricultural machinery power, crop sown area, number of agricultural employees, rural per capita agricultural output value and affected area of farmland. A regression analysis ensued. Following dimensionless processing of the original data, KMO and Bartlett tests were performed on the eight indicators. Results showed a KMO value of 0.690, meeting the factor analysis standard, and a P value of 0.000 (P < 0.05 indicating correlation between indicators); the findings rejected the null hypothesis, indicating the correlation between the indicators and the suitability of factor analysis.

Subsequently, principal component analysis of the index was conducted with the criterion that the eigenvalue exceeded 1 and the cumulative variance contribution rate was no less than 70%. The extracted components served as common factors, and the factor load matrix was further rotated through orthogonal rotation of maximum variance to obtain the factor score coefficient state.

According to Table 6, the first two variables from the principal component analysis summarize 93.2% of all independent variable information, demonstrating strong representativeness. The final principal component factor variables, denoted as FAC ₁ and FAC ₂, are obtained by performing maximum orthogonal variance rotation on these two variables. Table 6 shows the process and results of principal component analysis. Table 7 shows the mathematical relationship between principal component variables and initial variables. Then, the two principal component factors and the multiple regression analysis were carried out. The regression results revealed that the R ² value of the model was 0.943, coefficient P < 0.05, the model is significant indicating that the equation was well-fitted. The regression equation between lnCE and FAC ₁ and FAC ₂ can be obtained as follows:

(25)$$\ln CE = {-}0.071FAC_1 + 0.071FAC_2-1.252$$

Table 6. Conducting principal component analysis on influencing factors provided results on the contribution rates of factor variance

Ingredient represents the new variable extracted through principal component analysis, and characteristic root represents the variance of the principal component factor, representing the magnitude of the variance of the data interpreted by the principal component.

Table 7. Principal component score coefficient matrix of FAC ₁ and FAC ₂ shows the corresponding relationship between principal component factors and influencing factors

FAC ₁ and FAC ₂ are the final principal component factor variables; P, total agricultural output value (AOV) (100 million RMB); A, proportion of the agricultural output value among the total output value of agriculture, forestry, husbandry and fishery (AP) (%); T, urbanization rate, which is expressed as the proportion of the urban population among the total population of the BTH region (UR) (%); V, technical factor (i.e. the development level of agricultural technology, represented by the total power of agricultural machinery (AMP) (10⁴ kW)); U, is the sown area of crops (ACR) (hm²); C, number of agricultural employees (AE) (10 000 people); R, development level of the agricultural economy, expressed by the per capita agricultural output value of rural areas (ARV) (RMB/person); Y, the affected area of farmland (AFA) (hm²).

Equation (25) is the final regression equation obtained through principal component analysis. However, in order to more intuitively and obviously see the impact of influencing factors on carbon sinks, we replaced the principal component variable with the initial variable according to the mathematical relationship in Table 7, and obtained the regression equation between the initial variable and the carbon sink as follows:

(26)$$\eqalign{\ln ( CE) = & -1.252 + 0.022 \times \ln ( P) -0.002 \times \ln ( A) + 0.018 \cr & \quad \times \ln ( T) -0.009 \times \ln ( V) -0.013 \times \ln ( U) + 0.017 \cr & \quad \times \ln ( C) + 0.023 \times \ln ( R) -0.014 \times \ln ( Y) } $$

The obtained STIRPAT extended model of the net carbon sink in the BTH region is as follows:

(27)$$CE = 0.286 \times P^{0.022}A^{{-}0.002}T^{0.018}V^{{-}0.009}U^{{-}0.013}C^{0.017}R^{0.023}Y^{{-}0.014}$$

All letters in formulae (26–27) have the same meaning as in Eqn (24).

Robustness of the model refers to whether the evaluation and indicators can maintain the stable interpretation of the model when some parameters of the model are changed. Therefore, employing the Trimming method to assess the robustness of the STIRPAT model involves removing portions from both the beginning and end of the modelling data, shortening the time range, and observing whether the model exhibits any noticeable variations. The comparison of model parameters is shown in Table 8.

Table 8. Comparison between the Trimming method and the original model

P, total agricultural output value (AOV) (100 million RMB); A, proportion of the agricultural output value among the total output value of agriculture, forestry, husbandry and fishery (AP) (%); T, urbanization rate, which is expressed as the proportion of the urban population among the total population of the BTH region (UR) (%); V, technical factor (i.e. the development level of agricultural technology, represented by the total power of agricultural machinery (AMP) (10⁴ kW)); U, is the sown area of crops (ACR) (hm²); C, number of agricultural employees (AE) (10 000 people); R, development level of the agricultural economy, expressed by the per capita agricultural output value of rural areas (ARV) (RMB/person); Y, the affected area of farmland (AFA) (hm²); e, error term.

Upon review, it is evident that, compared to the original model, the symbols and parameters of the model after tail reduction treatment have shown no significant changes. The order of the degree of influence has also remained unchanged, with a significance level of 0.000, indicating a high level of model credibility.

Results

Temporal and spatial characteristics of net carbon sinks

Temporal and spatial characteristics of the net carbon sinks

Applying the method outlined in the previous section, we calculated the agricultural carbon sinks, carbon emissions and net carbon sinks for the 13 regions in the BTH region from 2009 to 2020. The changes in net carbon sinks during the study period are reported in Table 9.

Table 9. Carbon sink/carbon source change in Beijing-Tianjin-Hebei Region, 2009–2020

Statistical results in Table 9 and Fig. 2 show that the net carbon sink in the BTH region was 24 648 500 t in 2009, increasing to 34 592 900 t in 2020. The growth rate was 40.34%, with an average annual growth rate of 3.67%. Sequential growth was generally positive, except for slight decreases in 2014 and 2018. Agricultural carbon emissions underwent two stages: a stable period from 2009 to 2015, followed by a gradual decrease from 2016 onwards at a rate of 3.84%. Figure 3 illustrates the composition of agricultural carbon emissions, indicating an inflection point in 2016. Carbon emissions from agricultural land use, accounting for 53.37% of total emissions, showed high levels for nearly 10 years after 2009, decreasing annually post-2016. Livestock and poultry farming, the second-largest agricultural carbon source, accounted for 25.01% of emissions. Despite fluctuations, the overall trend was a slow decline. Crop respiration and straw burning, contributing smaller proportions, followed similar trends to livestock and poultry farming. Overall, increasing crop yield and continuous reduction in carbon emissions led to a rising agricultural carbon sink and net carbon sink.

Figure 2. Change of urban net carbon sink in Beijing-Tianjin-Hebei, 2009–2020.

Figure 3. Carbon emission composition of the Beijing-Tianjin-Hebei region from 2009 to 2020.

Analysing the change in the net carbon sink for each region, Shijiazhuang, Handan, Baoding, Hengshui and others contributed the most, while Qinhuangdao, Zhangjiakou, Chengde and others had the lowest net carbon sink. Spatially, high carbon sinks were mainly concentrated in the southern region, with generally low net carbon sinks in the north. Southern regions exhibited significantly higher crop yields than the north, and the eastern coastal region surpassed the west. Climate differences between the northwest and southeast regions, with the former being higher and the latter lower, contributed to this spatial variation. The North China Plain, located in the southeast and blocked by the Yanshan Mountains in the north, is not suitable for large-scale crop cultivation, resulting in significantly higher crop yields in the southern region. From the perspective of change, Tangshan, Langfang, Cangzhou and others maintained a stable net carbon sink, while all other regions, excluding Beijing, experienced varying degrees of increase.

Spatio-temporal variation of the net carbon sink intensity

The efficiency of the regional agricultural net carbon sink is reflected in the intensity of the net carbon sink. Illustrated in Fig. 4, the net carbon sink intensity in the BTH region rose from 261.01 kg/hm² in 2009 to 403.05 kg/hm² in 2020, displaying a fluctuating increasing trend. Apart from the years 2011 and 2017, the increasing trend remained relatively stable. This indicates continuous enhancement in the efficiency of the agricultural net carbon sink in the BTH region, with a significant increase in net carbon sink per unit of cultivated land. The net carbon sink intensity in the BTH region was categorized into five levels at four time points in 2009, 2012, 2016 and 2020.

Figure 4. Changes of net carbon sink intensity and carbon sink level in the Beijing-Tianjin-Hebei region from 2009 to 2020.

Analysing the change in Fig. 5, from a distribution perspective, the net carbon sink intensity in the southern region generally exceeded that in the northern region, indicating a tendency of aggregation. Based on net carbon sink intensity, the BTH region can be divided into three areas: (1) the northern region with low net carbon sink intensity, primarily including Zhangjiakou, Chengde, Qinhuangdao, etc.; (2) the region with medium net carbon sink intensity, mainly comprising coastal regions around the Bohai Rim, such as Tianjin, Tangshan, Cangzhou, etc.; (3) the region with high net carbon sink intensity in the south, encompassing Shijiazhuang, Handan, Xingtai and others. From the perspective of local region changes, the net carbon sink intensity of most regions increased, indicating significant overall improvement in the net carbon sink efficiency of the BTH region. However, the Chengde, Zhangjiakou and Qinhuangdao regions remained stable without a significant increase. Additionally, the net carbon sink intensity of Beijing exhibited a noticeable decreasing trend during the study period; it decreased from 257.40 kg/hm² in 2009 to 175.89 kg/hm² in 2020, marking a decline of 31.67%. This may be attributed to the direction of economic development. Agricultural development is not the primary focus of Beijing's economic direction; its crop output has experienced a significant decreasing trend, and most of its agricultural forms are small and micro-agriculture. This indicates a gradual decline in the scale and efficiency of agricultural production.

Figure 5. Change of urban net carbon sink intensity in Beijing-Tianjin-Hebei from 2009 to 2020. (a) 2009; (b) 2012; (c) 2016; (d) 2020.

Spatiotemporal variation of carbon sink levels

The agricultural carbon sink level reflects the regional structure of the agricultural carbon sink/source, with a higher ratio indicating a higher carbon sink level. As indicated in Fig. 4, the agricultural carbon sink level in the BTH region has shown a consistent and increasing trend, signifying a greater change range in the carbon sink compared to carbon emissions. The accounting results of the net carbon sink reveal a decreasing trend in agricultural carbon emissions in the BTH region, coupled with a continuous increase in the carbon sink due to improved crop yield. This forms the fundamental reason for the ongoing enhancement of the carbon sink level.

Examining the 12-year change in carbon sink level in the BTH region, as illustrated in Fig. 6, the alterations in agricultural carbon sink level and net carbon sink intensity align closely. The southern region exhibits a greater degree of change in carbon sink level compared to the northern region. Handan and Cangzhou display the most noticeable changes, experiencing a significant increase in their carbon sink levels during the study period, although these levels significantly differ from those observed in the southern region. Langfang undergoes minimal change, while Beijing demonstrates a decreasing trend. Over time, the BTH region forms a circular distribution with Beijing as the centre, showcasing an increase in carbon sink level from the north to the south. This pattern suggests well-developed agricultural carbon sinks in the southern region, while the northern region maintains a stable state. The ratio of carbon sink to carbon emission in Beijing declined from 2.17 in 2009 to 1.76 in 2020, making it the sole region in the BTH region with a decreasing carbon sink level. In general, regional disparities are the primary drivers of changes and variations. Therefore, developing agriculture with local characteristics according to regional conditions has become a crucial strategy to address this situation. Due to its unique urban positioning, Beijing's agricultural production is gradually diminishing. In this scenario, it should compensate for the ‘carbon sink’ from alternative perspectives to mitigate the potential risk of a ‘carbon deficit’.

Figure 6. Change of urban carbon sink level in Beijing-Tianjin-Hebei from 2009 to 2020. (a) 2009; (b) 2012; (c) 2016; (d) 2020.

Spatial autocorrelation analysis based on Moran's index

Global spatial autocorrelation analysis

The Moran's I Index of the net carbon sink in the BTH region during 2009–2020 was calculated, and the results are reported in Table 10.

Table 10. The spatial autocorrelation analysis results of agricultural net carbon sink in the BTH region

Moran's I value is the slope of the line that best fits the relationship between neighbouring income values and each polygon's income in the dataset; Z score describes how far away from the mean (or average) the data lies in a normally distributed sample.

As indicated in Table 10, the Moran's I Index for the total net carbon sink in the BTH region consistently showed positive values, z is a multiple of the standard deviation, accompanied by P values below 0.05 (indicating a successful 95% confidence test). This suggests the presence of autocorrelation among the net carbon sinks from 2009 to 2020, signifying a prevalent trend of mutual influence and interaction within the net carbon sinks in the BTH region. The Moran's I Index values ranged between 0.393 and 0.542, averaging 0.476. Despite fluctuations in the Moran's I Index during the latter part of the study period, an overall increasing trend persisted, reaching its peak at 0.581 in 2017, showcasing the strongest spatial autocorrelation. This trend illustrates an ongoing enhancement in the spatial aggregation of net carbon sinks in the BTH region over time.

Local spatial autocorrelation analysis

The Local Moran's Index was employed to depict the spatial clustering of net carbon sinks in the BTH region, with Local Indicators of Spatial Autocorrelation (LISA) cluster maps generated for the years 2009, 2012, 2016 and 2020 (refer to Fig. 7). The analysis reveals a distinct local spatial aggregation phenomenon and a noticeable differentiation pattern in the net carbon sinks across the BTH region. The highly aggregated areas are predominantly situated in the southern region, prominently featuring Shijiazhuang, Hengshui and Xingtai. These regions, characterized by a warm climate, serve as crucial agricultural production hubs in the BTH region, showcasing substantial advantages in agricultural industry development. Conversely, Chengde stands as the core of the L-L cluster area in the north, where the climate tends to be colder, and agricultural production conditions are relatively weaker. In terms of development and change, the highly aggregated area assumes a radiating and driving role, with its clustering remaining consistent over the study period and exhibiting a weak radiating driving effect on surrounding regions. In contrast, the low clustering area, centred around Chengde, undergoes continuous changes throughout the study period, expanding to surrounding areas. Tangshan only displayed H-L clustering in 2009, and over time, the agglomeration effect weakened, possibly influenced by surrounding regions and exhibiting a regional development-oriented industrial structure. Hebei, characterized by diverse landforms and topographic features, has scattered mountains, grasslands and hills in the northwest and north, lacking spatial uniformity. However, the southern part of Hebei, predominantly plains, is suitable for large-scale agricultural development, displaying better spatial integrity and a relatively concentrated spatial aggregation effect than the northern part.

Figure 7. Local Indicators of Spatial Autocorrelation (LISA) cluster map of urban net carbon sinks in the Beijing-Tianjin-Hebei (BTH) region, 2009–2020. (a) 2009; (b) 2012; (c) 2016; (d) 2020.

Analysis of inequality

Per capita carbon sink Theil index regional inequality

We calculated the Thiel indexes T(G) and T(P) with agricultural GDP and agricultural population as weights. These indexes reflect the correlation between the economic development level and agricultural population size with the agricultural net carbon sink. The study area was divided into five regions for analysing regional differences: eastern Hebei, northern Hebei, central Hebei, southern Hebei and Beijing-Tianjin.

In Table 11, the T(P) index consistently remained above 0.102, with an average of 0.122, exhibiting a gradual ‘N’-shaped change. This suggests that the per capita net carbon sink in the BTH region showed significant regional differences with fluctuating trends from 2009 to 2020. According to Figs 8(a) and (d), the ratio of intra- and inter-regional difference contribution rates for the per capita net carbon sink was 15.7: 84.3% in 2009 and 27.2: 72.8% in 2020. This indicates that the overall difference in the per capita net carbon sink in the BTH region was mainly due to inter-regional differences. Among intra-regional differences, the five regions had varying impacts. The central Hebei region had the highest contribution rate, averaging 29.5%, while Beijing, Tianjin and southern Hebei contributed 22.1 and 18%, respectively. The regional differences showed a steady increasing trend, except for a significant fluctuation in 2010.

Table 11. Changes in Theil index from 2009 to 2020

T(P), per capita carbon sink Theil index; T(G), the carbon sink intensity Theil index.

Figure 8. Change of T(P) and T(G) contribution rates in BTH cities from 2009 to 2020. (a) T(P) contribution rate of intra-regional and inter-regional differences; (b) T(P) regional contribution rate; (c) T(G) contribution rate of intra-regional and inter-regional differences; (d) T(G) regional contribution rate.

Carbon sink intensity Theil index of regional inequality

Table 11 shows that the T(G) index consistently remained above 0.07 throughout the study period, with an average of 0.083. The peak occurred in 2009, followed by a linear decreasing trend. By 2020, the index had decreased to 0.070 with noticeable regional differences. This decreasing trend was more pronounced than that of the T(P) index. When examining Figs 8(c) and (d) and decomposing regional differences, the ratio of intra- and inter-regional difference contribution rates for the net carbon sink intensity of agricultural output value was 15.7: 84.3% in 2009 and changed to 77.1: 22.9% in 2020. This suggests a shift from intra-regional to inter-regional differences in the net carbon sink intensity of the agricultural output value in the BTH region. The intra-regional contribution rate of each region showed continuous increase. Central Hebei had the highest contribution rate, averaging 23.0%, while southern and eastern Hebei contributed 11.9 and 3.79%, respectively. The contribution rate of regional differences increased, with the most noticeable change observed in Beijing, Tianjin and southern Hebei.

Overall, the mean values of T(P) and T(G) over the 12-year study period were 0.122 and 0.070, respectively. This indicates that the Thiel index T(P) weighted by agricultural population was larger than T(G) weighted by agricultural GDP. T(P) better reflects regional differences in the agricultural net carbon sink in the BTH region. However, the matching degree between regional net carbon sink and agricultural output value was higher than that between agricultural population size, and agricultural population is a significant factor contributing to regional differences in net carbon sinks. As can be seen from Table 11, over the past 12 years, T(P) and T(G) values showed an overall declining trend. This suggests that the severe imbalance in agricultural development in Beijing, Tianjin and Hebei has improved, but there are still noticeable coordination challenges among different regions that need further improvement.

Driver analysis based on the extended STIRPAT model

According to the regression equation of Eqn (27), the degree of influence of each factor on the net carbon sink in the BTH region was as follows: ARV > AOV > UR > AE > AFA > ACR > ATL > AP. Each 1% change in AOV, AP, UR, AMP, ACR, AE, ARV and AFA will result in a change of 0.022, 0.002, 0.018, 0.009, 0.013, 0.017, 0.023 and 0.014% in net carbon sink. Among them, AOV, UR, AE and ARV were found to have significant promoting effects on the net carbon sink, among which the promoting effect of ARV was the strongest. AP, AMP, ACR and AFA had reverse inhibitory effect on agricultural net carbon sink, among which the inhibiting effect of AFA area was the strongest. AP was the least restrictive.

The total agricultural output value and rural per capita output value were found to be the main contributors to the increase of the agricultural net carbon sink, which indicates that the rural economic development level is the main factor that promotes the agricultural net carbon sink. The total agricultural output value and per capita rural output value in the BTH region continued to grow from 2009 to 2020. This indicates that the continuous growth of the rural economic level led to the increase of the agricultural industry input and the promotion of the development of agriculture towards low-carbon agriculture. The urbanization level and the number of agricultural employees were found to be secondary factors that promote the net carbon sink. The urbanization level of BTH increased from 50.19% in 2009 to 64.83% in 2020, indicating that urban expansion does not inhibit the net carbon sink effect, but instead promotes the net carbon sink effect. In addition, affected by the increase of the urbanization rate, the proportion of the rural population decreased continuously during the study period, and the number of agricultural employees in the BTH region presented an obvious declining trend. According to the analysis results, the number of agricultural employees is a positive factor that promotes the growth of the net carbon sink. Therefore, to compensate for the negative effect brought by the decline of agricultural employees, agricultural professionals should be vigorously developed to achieve the efficient management of agricultural production to offset the consequences of urbanization.

The sown area of crops and the affected farmland area are key factors hindering the agricultural net carbon sink. In the BTH region, the crop sown area exhibited a declining trend, contrary to the net carbon sink trend. The increase in net carbon sink intensity compensates for the negative impact of reduced crop sown area. The relationship between affected farmland area and net carbon sink aligns with expectations, as farmland disasters directly impact crop yield and, subsequently, the carbon sink. The rural technical level is a significant inhibiting factor. Improved agricultural mechanization has not positively contributed to increased carbon sink effectiveness, primarily due to the overuse of machinery and the resulting energy consumption, leading to higher carbon emission intensity. Therefore, it is crucial to address excessive carbon emissions resulting from agricultural mechanization.

Discussion

Issues facing agricultural development and carbon neutrality

Based on the aforementioned research findings, we can summarize the key trends in the development of agricultural carbon sinks in the BTH region. In recent years, there has been notable improvement in both the capacity and effectiveness of carbon sinks, significantly contributing to the advancement of agricultural development (Sun et al., 2023). However, it is crucial to acknowledge and address the challenges accompanying this positive trend. The predominant challenge lies in the current regional development imbalance. The agricultural carbon pooling in the BTH region exhibits a discernible pattern of north-south differentiation and aggregation. Over time, this disparity has shown a tendency to widen, closely linked to the overall imbalance in agricultural development (Gan et al., 2023). According to the research conducted by Kong and Cheng (2017), Hebei surpasses Beijing and Tianjin in both land output rate and its growth rate. Additionally, the primary industry's proportion in Hebei is considerably higher than that in Beijing and Tianjin, highlighting substantial developmental differences. Nevertheless, within Hebei, internal disparities are also pronounced. Li et al. (2019) discovered, through their study on agricultural green factor productivity in Hebei, that hot spots have predominantly expanded to central and eastern Hebei, while cold spots have shifted from northern Hebei to the two wings. This spatial concentration and distribution of hot and cold spots persist over time. Excessive industrial concentration poses challenges to the coordinated development of these regions. Various factors contribute to this phenomenon. Scholars argue that the primary reasons behind the agricultural coordinated development challenges in Beijing, Tianjin and Hebei include the inequality in regional administrative status, inconsistency in development interests and the absence of an effective cooperation mechanism (Guo et al., 2017). Given Beijing's role as the political and cultural centre at the core of the BTH region, the coordinated development plan primarily positions Tianjin and Hebei to support the capital's economic development on the supply side (Fan et al., 2022), Consequently, the agricultural functions of the central and southern Hebei regions become more pronounced (Jiao et al., 2020). Despite China's efforts to promote coordinated agricultural development in the BTH region, there is still ample room for improvement (Xiao et al., 2022).

Provide guidance based on the challenges faced

In order to ensure the efficient and stable development of carbon sink capacity in the BTH region and realize regional collaborative integration as soon as possible, the following development suggestions are put forward according to the development characteristics and functional positioning of the BTH region and combined with the research results:

1. (1) According to the regional agricultural production characteristics, agricultural emissions reduction and industrial structure upgrading should be carried out according to the local conditions.

The government should pay attention to the unbalanced development of regional agricultural production. For areas with high net carbon sinks, focus should be placed on carbon emission reduction, transition to low-carbon agriculture, industrial layout planning and the reduction of agricultural inputs such as pesticides and fertilizers. In Hebei, the proportion of straw burning should be reduced, and carbon emission targets should be constrained. In areas with a low net carbon sink, the industrial structure should be optimized, the agricultural production cost input should be increased, the agricultural production efficiency should be improved and characteristic and ecological agriculture should be developed according to the local conditions. Moreover, the agricultural industrial structure should be avoided, the further widening of regional differences should be prevented and the collaborative integration of high-carbon-sink demonstration areas should be jointly established by Beijing, Tianjin and Hebei.

1. (2) The construction of agricultural regions with a high carbon sink and regional coordinated development should be promoted.

High-carbon-sink areas should play a leading and exemplary role for the radiating surrounding areas, thus driving the development of regional rural low-carbon industries and realizing the balance of regional agricultural industries. For regions with a declining agricultural carbon sink capacity, the carbon sink/source structure should be adjusted according to the industry type to avoid a carbon deficit. Via the driving effect of the high-carbon-sink regions in the BTH region on the surrounding areas, support can be provided for the integration of agriculture in the region.

1. (3) The construction of new energy in rural areas should be promoted and the energy efficiency should be improved.

The technical level of rural areas has a negative effect on the net carbon sink, indicating that mechanized agricultural production leads to a serious carbon emission problem. It is necessary to strengthen technological innovation and industrial upgrading, innovate the technical level of agricultural mechanization and promote the development and application of energy-saving and environmental protection technologies, e.g. the use of clean energy such as electric energy and solar energy. Moreover, the proportion of fossil fuels such as diesel should be reduced, and importance should be attached to the efficiency of agricultural resource utilization. The waste of ineffective resources should also be reduced.

1. (4) The construction of high-quality agricultural talents should be promoted and resource management should be coordinated.

The government should pay attention to the negative impact of the improvement of the urbanization level on the decrease of the number of rural employees. It should also promote the construction of new urbanization, the centralized management of the agricultural production process, the development of high-quality and professional agricultural talents and the centralized planning of rural surplus labour. Moreover, it should aggregate the management of the population, land and other production resource factors to achieve efficient agricultural production. This will provide a solid guarantee for the construction of the BTH high-carbon sink demonstration zone, as well as strong support for the regional agricultural integration of BTH by establishing an agricultural demonstration zone for population-intensive industries.

Most of the data sources used in this research method are from statistical yearbooks, and the statistical types are relatively limited. For counties and smaller research scales, the proportion of characteristic agricultural industry is higher than that of main agriculture in some regions, and the statistical data of statistical yearbooks cannot reflect the degree of regional agricultural development. Therefore, this method is suitable for research at the provincial and municipal scale. It should be optimized according to the local agricultural development conditions. In addition, most of the parameters used in this study are summarized from previous research results and have universal applicability. However, for research on a small regional scale, these parameters may be different from the actual local situation. If detailed estimates are required, these parameters should be adjusted according to the actual situation in the study area.

Limitations and future prospects

In addressing carbon emissions within agricultural production, this study takes into account the dual attributes of carbon sink/carbon source. It avoids the singular perspective of carbon storage or emissions found in previous research and calculates the effects of carbon sink, carbon source and net carbon sink in agricultural production from multiple angles. The study analyses the carbon sink effect and its change trends in China's typical agricultural production agglomeration areas over time and space, aiming to provide guidance for the development of regional green agriculture with low-carbon production practices.

It is crucial to note that the focus of this study on carbon sink and carbon emission pertains to the storage and release of carbon in the agricultural production process. However, it is well-known that carbon undergoes constant recycling in nature, and most of the carbon stored in crops will eventually return to nature through biological consumption. Therefore, preserving carbon reserves as much as possible during the consumption process remains an important issue we must confront.

In the realm of agricultural production, reducing carbon emissions and increasing the ratio of carbon sink to carbon source are also critical considerations. Fortunately, our research indicates that crop yields and carbon sinks continue to increase in agricultural activities in the BTH region between 2009 and 2020, while carbon emissions show a downward trend. This suggests that the carbon sink effect of agricultural production is moving in a positive direction, offering valuable insights for the establishment of low-carbon agriculture.

Conclusions

The net carbon sink in the BTH region demonstrates a gradual increasing trend, with a concurrent decline observed in agricultural carbon emissions. Agricultural land use emissions and livestock farming emissions stand out as the primary sources of carbon emissions. Meanwhile, the efficiency of agricultural net carbon sequestration and the level of agricultural carbon sink are continually improving. The regional net carbon sink exhibits evident spatial autocorrelation, with notable north–south disparities and pronounced spatial inequalities. The level of rural economic development emerges as a key factor promoting the increase in agricultural net carbon sequestration. In summary, from 2009 to 2020, both the quantity and level of carbon sequestration in the BTH region have significantly increased. However, regional development imbalances persist, emphasizing the importance of coordinated and integrated regional development as a crucial direction for achieving agricultural carbon neutrality.