Ecological system carbon reserve remote sensing calculation method considering earth surface type and quality

Abstract

The invention provides an ecological system carbon reserve remote sensing measurement and calculation method considering earth surface type and quality, which comprises the following steps of: importing an average carbon density table, a maximum carbon density table and a minimum carbon density table of a research area and land utilization data of the research area into a carbon storage and storage module of InVEST software; a theoretical original carbon reserve grid image # imgabs0 #, a theoretical maximum original carbon reserve grid image # imgabs1 # and a theoretical minimum original carbon reserve grid image # imgabs2 # of each land class of a research area are calculated through a built-in algorithm; a formula 30 is used for calculating a kth part carbon reserve grid image # imgabs3 # in a certain land class i; a formula 31 is used for calculating the total carbon reserve of the land class i of the research area; ai represents the area of the land class i; the sum of the total carbon reserves of the districts is the sum of the total carbon reserves of the research area.

Description

Technical Field

[0001] The present application relates to remote sensing technology, and in particular to a remote sensing calculation method for ecosystem carbon storage that takes into account both surface type and quality. Background Art

[0002] Greenhouse effect and climate warming are global ecological problems faced by mankind. Studies have shown that a 2°C increase in global average temperature will lead to a series of irreversible hazards such as frequent extreme climate, reduced food production, and species extinction, which will seriously threaten human survival and development. Wayne S et al.'s research shows that the top 25 contributors to the unrealized additional land carbon storage potential account for nearly three-quarters of the global total, with China ranking fourth with a potential of more than 15Pg C. Effectively exerting the carbon sequestration role of forests, grasslands, wetlands, oceans, soils and permafrost to increase the increment of ecosystem carbon sinks is a powerful measure. As the concept of ecological restoration and protection gradually takes root in the hearts of the people, a series of laws and regulations have been formulated to repair the ecological environment. Carrying out ecological construction work such as land remediation and planting of forest and grass vegetation in areas damaged by human activities will inevitably significantly improve vegetation coverage and change the carbon sequestration capacity of the ecosystem. Therefore, scientifically and accurately assessing the current status of carbon storage and exploring its spatial distribution and the influence of driving factors have a certain guiding role in the formulation of ecological protection strategies and the implementation of carbon sink enhancement measures. Current research on carbon storage measurement mainly focuses on other terrestrial ecosystems such as forests and grasslands, as well as urban ecosystems that are more closely related to human survival and development. In the long run, it is necessary to evaluate the effectiveness of ecological governance from the perspective of carbon storage to provide a reference for future ecological governance and decision-making.

[0003] As the global climate change problem becomes increasingly serious, there have been multi-faceted and multi-angle studies focusing on the measurement and assessment of carbon storage. Foreign research on carbon storage has been going on for about 60 years. Since the 1950s, international scholars have conducted research on the estimation of global soil carbon storage. Foreign research on carbon storage mainly includes carbon sources, carbon sinks and carbon storage measurement; research methods mainly include sample plot inventory method, bookkeeping model, InVEST (Integrated Valuation of Ecosystem Services and Trade-offs) model, etc.; the research scale is mainly focused on the impact of land use changes on the overall carbon cycle over a long period of time. my country's research on carbon storage started later than that of foreign countries. Since the 1990s, Chinese scholars have calculated the carbon storage of soil and vegetation based on soil census and forest resource survey data, and further analyzed the role of land use changes in the carbon cycle of ecosystems. In China, research content mainly focuses on carbon storage measurement and exploring the impact of land use changes on carbon storage; research methods also include sample site inventory, process models and remote sensing models, among which the InVEST model is used to measure carbon storage. In the selection of research areas, more attention is paid to functional areas such as forests, grasslands, wetlands, and ecologically fragile areas such as arid areas and karst areas.

[0004] Regarding the specific research hotspot of carbon storage measurement, the following problems were found in the current research: InVEST model is widely used in carbon storage measurement using models. Researchers often use carbon density in the form of a table to participate in the process of calculating carbon storage using the InVEST model, that is, the carbon density of different land use types is different, and the carbon density of the same land use type is consistent. However, the actual situation is that the same type of land has differences in quality, and its carbon density should not be regarded as equal. At this time, the impact of surface quality on carbon density should be considered, so that it can reflect the land use type and present the differences in quality of the same surface type in the process of carbon storage measurement. Summary of the invention

[0005] In view of the above problems, this application aims to propose a remote sensing calculation method for ecosystem carbon storage that takes into account both surface type and quality. In the process of carbon storage measurement, it reflects the land use type while presenting the differences in quality of the same surface type.

[0006] The remote sensing calculation method of ecosystem carbon storage taking into account the surface type and quality of the present application includes:

[0007] The average carbon density table, maximum carbon density table, minimum carbon density table and land use data of the study area were imported into the carbon storage and sequestration module of InVEST software, and the theoretical original carbon storage raster image of each type of land in the study area was calculated through its built-in algorithm. Theoretical maximum original carbon storage raster image and theoretical minimum original carbon storage raster images The resolution of the imagery is consistent with the land use data;

[0008] use To calculate the carbon storage grid image of the kth part in a certain land type i

[0009] Among them, k represents a part of carbon storage, k represents the total carbon storage when it is total, k represents the aboveground carbon storage when it is above, k represents the underground carbon storage when it is below, k represents the soil carbon storage when it is soil, and k represents the dead organic carbon storage when it is dead; i represents the number of each type of land, 1 for cultivated land, 2 for forest land, 3 for grassland, 4 for water area, 5 for building land, and 6 for unused land; CDI is the carbon density index, which is used to indicate the surface quality;

[0010] is the aboveground carbon stock of land type i calculated using formula 30; is the underground carbon stock of land type i calculated using formula 30; is the carbon stock in the soil of land type i calculated using formula 30; is the dead organic carbon stock of land type i calculated using formula 30;

[0011] Then use To calculate the total carbon storage of land type i in the study area, A i Represents the area of ​​land type i; the sum of the total carbon storage of each land type is the total carbon storage of the study area.

[0012] Preferably, CDI is based on Calculation is performed; kNDVI is the greenness index; kSVMI is the humidity index;

[0013] in,

[0014]

[0015] Among them, tanh is the hyperbolic tangent function, sgn is the sign function, and ρ NIR is the near-infrared band, ρ Red It is the red light band;

[0016]

[0017] In the formula, ρ Blue is the blue light band, ρ Red For the red light band, is the short-wave infrared band, and j can be 1, 2 and 3, representing short-wave infrared bands of different wavelength ranges.

[0018] Preferably,

[0019]

[0020] in, is the soil carbon density obtained based on the annual maximum precipitation HAP, and are the biomass carbon density obtained according to the annual maximum temperature HAT and the annual maximum precipitation HAP; C′ i and C″ i Data represent the study area and the whole country, respectively; is the maximum biomass carbon density correction factor; is the maximum soil carbon density correction factor

[0021] Use the maximum biomass carbon density correction factor Multiply them with the aboveground carbon density, underground carbon density and dead organic matter carbon density columns in the national carbon density data, and use the maximum soil carbon density correction factor Multiply it with the soil carbon density column in the national carbon density data to obtain the maximum carbon density table of the study area.

[0022] Preferably,

[0023]

[0024]

[0025] in, is the soil carbon density obtained based on the annual minimum precipitation LAP, and are the biomass carbon density obtained according to the annual minimum temperature LAT and the annual minimum precipitation LAP; C′ i and C″ i Data represent the study area and the whole country, respectively; is the minimum biomass carbon density correction factor, is the minimum soil carbon density correction factor;

[0026] Use minimum biomass carbon density correction factor Multiply them with the aboveground carbon density, underground carbon density and dead organic matter carbon density columns in the national carbon density data, and use the minimum soil carbon density correction factor Multiplying it with the soil carbon density column in the national carbon density data, the minimum carbon density table of the study area was obtained.

[0027] Preferably,

[0028] C SP =3.3968×MAP+3996.1(R 2 =0.11) (Formula 2);

[0029] C BP =6.798×e 0.0054×MAP (R 2 =0.70) (Formula 3);

[0030] C BT =28×MAT+398(R 2 =0.47, P<0.01) (Formula 4);

[0031]

[0032] K B =K BT ×K BP (Formula 7);

[0033]

[0034] Among them, C SP is the soil carbon density obtained based on the average annual precipitation MAP, C BT and C BP are the biomass carbon density obtained based on the average annual temperature MAT and the average annual precipitation MAP;

[0035] C i ′ and C i ′′ represents the data of the study area and the whole country respectively. When it means the biomass carbon density is improved according to precipitation, i=BP; when it means the biomass carbon density is improved according to temperature, i=BT; when it means the soil carbon density is improved according to precipitation, i=SP; K B is the biomass carbon density correction factor; K S is the soil carbon density correction factor; the biomass carbon density correction factor K is used B Multiply them with the aboveground carbon density, underground carbon density and dead organic matter carbon density columns in the national carbon density data, and use the soil carbon density correction factor K S Multiplying it with the soil carbon density column in the national carbon density data, we get the average carbon density table for the study area.

[0036] The remote sensing calculation method for ecosystem carbon storage in the present application takes into account both the surface type and quality. In the process of carbon storage calculation, it reflects the land use type while presenting the differences in quality of the same surface type, making the carbon storage calculation more accurate. BRIEF DESCRIPTION OF THE DRAWINGS

[0037] Figure 1 This is a flow chart of the remote sensing calculation method for ecosystem carbon storage that takes into account surface type and quality in this application.

[0038] Figure 2 This is an overview map of the Ningdong mining area.

[0039] Figure 3 This is a graph showing the results of carbon reserve estimation in the mining area based on the InVSET model.

[0040] Figure 4 This is a graph showing the results of carbon reserve estimation in mining areas based on the InVSET model, taking both type and quality into consideration.

[0041] Figure 5 This is the spatial distribution map of the improved temporal correlation between aboveground carbon storage and MYD15A2H products from 2017 to 2022. DETAILED DESCRIPTION

[0042] The calculation of original carbon storage is achieved using the InVEST model. The InVEST model is a comprehensive assessment model of ecosystem services and trade-offs jointly developed by Stanford University, The Nature Conservancy and other institutions. The carbon storage and sequestration module can use land use to estimate the carbon storage in the current landscape. Carbon storage mainly includes four parts, namely aboveground biomass carbon pool, underground biomass carbon pool, dead organic matter carbon pool and soil carbon pool.

[0043] In this application, national carbon density data, land use data (Land Use and Land Cover, LULC) and climate data are used to calculate the theoretical maximum and minimum carbon storage corresponding to the study area. Then, the surface greenness index kNDVI (kernel Normalized Difference Vegetation Index) and the moisture index kSVMI (kernel-based optical remote sensing Index of Soil and Vegetation Moisture) are introduced, and the two are coupled by arithmetic average to generate the carbon density index CDI (Carbon Density Index) to linearize the original carbon storage calculated by the traditional method, and add the surface mass to the carbon storage calculation process.

[0044] 1. Calculation of carbon storage in original ecosystems

[0045] The calculation formula for the original ecosystem carbon storage is:

[0046]

[0047] In the formula, i represents the number of each land type, and cultivated land, forest land, grassland, water area, construction land and unused land are 1-6 respectively. i is the area of ​​a certain land type, in hectares (hm 2 ). Refers to the sum of carbon stocks of a certain land type calculated using the original method (OM). Refers to the aboveground biomass carbon pool of a certain land type. Refers to the underground biomass carbon pool of a certain land type. Refers to the dead organic carbon pool of a certain land type. Refers to the soil carbon pool of a certain land type, and the unit is ton (t). Taking cultivated land as an example, the total carbon storage of cultivated land calculated using the original method is The aboveground biomass carbon pool of cultivated land is The carbon pool of underground biomass in cultivated land is The dead organic matter carbon pool in cultivated land is The soil carbon pool of cultivated land is The area of ​​cultivated land in the study area is A1, and so on.

[0048] To obtain each carbon pool, we must first obtain a carbon pool raster image in pixels, which depends on land use data and carbon density data. General research usually corrects the national carbon density obtained from the survey according to local climate conditions in the study area, such as temperature and precipitation. The correction of aboveground biomass carbon density refers to the research of Chen Guangshui, as shown in Formula 4; the correction of soil carbon density and biomass carbon density refers to the research of Alam and Giardina, as shown in Formula 2 and Formula 3. The parameters involved include mean annual precipitation (MAP) and mean annual temperature (MAT):

[0049] C SP =3.3968×MAP+3996.1(R 2 =0.11) (Formula 2)

[0050] C BP =6.798×e 0.0054×MAP (R 2 =0.70) (Formula 3)

[0051] C BT=28×MAT+398(R 2 =0.47,P<0.01) (Formula 4)

[0052] In the formula, the average annual precipitation MAP (mm) comes from the ERA5-LandMonthlyAggregated-total_precipitation band in the ERA5 meteorological reanalysis data, and the mean annual temperature (MAT) comes from the ERA5-Land Monthly Aggregated-temperature_2m band in the ERA5 meteorological reanalysis data. The average of all monthly data is used to obtain MAP and MAT. SP is the soil carbon density (t / hm2) obtained based on the average annual precipitation MAP (mm) 2 ),C BT and C BP are the biomass carbon density (t / hm2) obtained based on the average annual temperature MAT (℃) and the average annual precipitation MAP (mm). 2 ).

[0053]

[0054] K B =K BT ×K BP (Formula 7)

[0055]

[0056] Where C′ i and C″ i Represent the data of the study area and the whole country respectively. When the biomass carbon density is improved according to precipitation, i = BP; when the biomass carbon density is improved according to temperature, i = BT; when the soil carbon density is improved according to precipitation, i = SP. According to formula 7 and formula 8, the biomass carbon density correction coefficient K can be obtained. B and soil carbon density correction factor K S The national carbon density data obtained from the survey is known, as shown in Table 1. The above formula is used to improve the carbon density according to the actual climate conditions of the study area. The specific operation process is: Use the biomass carbon density correction coefficient K B The aboveground carbon density, underground carbon density and dead organic matter carbon density in the study area were obtained by multiplying them with the aboveground carbon density, underground carbon density and dead organic matter carbon density in Table 1. The soil carbon density correction factor K was used. s The local soil carbon density of the study area was obtained by multiplying it with the soil carbon density column in Table 1, and finally the average carbon density table of the study area was obtained.

[0057] Table 1 National carbon density data obtained from the survey

[0058]

[0059] The annual average temperature and annual average precipitation in Formula 2-Formula 8 are replaced by the annual maximum precipitation (Highest Annual Precipitation, HAP) and the annual maximum temperature (Highest Annual Temperature, HAT), and the theoretical maximum carbon density of the study area after preliminary improvement is obtained. The formulas involved are shown in Formula 9-Formula 15:

[0060]

[0061] In the formula, the Highest Annual Precipitation (HAP) comes from the ERA5-Land Monthly Aggregated-total_precipitation band in the ERA5 meteorological reanalysis data, and the Highest Annual Temperature (HAT) comes from the ERA5-Land MonthlyAggregated-temperature_2m band in the ERA5 meteorological reanalysis data. The maximum value of all monthly data can be used to obtain HAP and HAT. is the soil carbon density (t / hm2) obtained based on the annual maximum precipitation HAP (mm) 2 ), and are the biomass carbon density (t / hm2) obtained based on the annual maximum temperature HAT (℃) and the annual maximum precipitation HAP (mm). 2 ). C′ i and C″ i Represent the data of the study area and the whole country respectively, as above. According to formula 14 and formula 15, the maximum biomass carbon density correction coefficient can be obtained and soil carbon density correction factor Use the maximum biomass carbon density correction factor Multiply them with the aboveground carbon density, underground carbon density and dead organic matter carbon density columns in the national carbon density data, and use the maximum soil carbon density correction factor Multiply it with the soil carbon density column in the national carbon density data to obtain the maximum carbon density table of the study area.

[0062] Replacing the annual average temperature and annual average precipitation in Formula 2-Formula 8 with the annual minimum precipitation (Lowest Annual Precipitation, LAP) and the annual minimum temperature (Lowest Annual Temperature, LAT), the preliminary improved theoretical minimum carbon density of the study area is obtained, as shown in Formula 16-Formula 22.

[0063]

[0064]

[0065] In the formula, the annual minimum precipitation (Lowest Annual Precipitation, LAP) comes from the ERA5-Land Monthly Aggregated-total_precipitation band in the ERA5 meteorological reanalysis data, and the annual minimum temperature (Lowest Annual Temperature, LAT) comes from the ERA5-Land MonthlyAggregated-temperature_2m band in the ERA5 meteorological reanalysis data. The minimum value of all monthly data is used to obtain LAP and LAT. is the soil carbon density (t / hm2) obtained based on the annual minimum precipitation LAP (mm) 2 ), and are the biomass carbon density (t / hm2) obtained based on the annual minimum temperature LAT (℃) and the annual minimum precipitation LAP (mm). 2 ). C′ i and C″ i Represent the data of the study area and the whole country respectively, as above. According to formula 21 and formula 22, the minimum biomass carbon density correction coefficient can be obtained and soil carbon density correction factor Use minimum biomass carbon density correction factor Multiply them with the aboveground carbon density, underground carbon density and dead organic matter carbon density columns in the national carbon density data, and use the minimum soil carbon density correction factor Multiplying it with the soil carbon density column in the national carbon density data, the minimum carbon density table of the study area was obtained.

[0066] Import the above-mentioned average carbon density table, maximum carbon density table, minimum carbon density table and land use data of the study area into the carbon storage and sequestration module of InVEST software, and then calculate the theoretical original carbon storage grid image of each type in the study area every year through the built-in algorithm of the software. Theoretical maximum original carbon storage raster image and theoretical minimum original carbon storage raster images The resolution of the image is consistent with the land use data. k represents a part of the carbon storage. When k is total, it represents the total carbon storage. When k is above, it represents the aboveground carbon storage. When k is below, it represents the underground carbon storage. When k is soil, it represents the soil carbon storage. When k is dead, it represents the dead organic carbon storage. i represents the number of each type of land. Cultivated land is 1, forest land is 2, grassland is 3, water area is 4, construction land is 5, and unused land is 6. Using the obtained raster image and formula 1, the carbon storage and total carbon storage of each part of the ecosystem in the study area can be calculated using the original method.

[0067] 2. Calculation of Carbon Density Index

[0068] The carbon density index CDI is a coupling of the greenness index kNDVI and the moisture index kSVMI, and is used to indicate the surface quality of the ecosystem in the study area.

[0069] (1) Greenness Index kNDVI

[0070] Greenness indices such as NDVI and EVI can reflect the differences in vegetation coverage and growth conditions to a certain extent, and are particularly suitable for monitoring vegetation with vigorous growth and high coverage. In 2021, Camps-Valls et al. proposed a nonlinear vegetation index kNDVI, which maximizes the use of spectral information and solves long-standing problems in satellite earth observations of the terrestrial biosphere. At the same time, it can more accurately measure terrestrial carbon sources, carbon sink dynamics, and the potential to stabilize atmospheric CO2 and mitigate global climate change. Under various environmental conditions, such as dense woodlands, grasslands, mixed forests, etc., kNDVI shows stronger stability and robustness than traditional indices such as NDVI, so it was selected. ini As a greenness index to improve carbon density. Its calculation formula is:

[0071]

[0072] Where tanh is the hyperbolic tangent function, sgn is the sign function, and ρ NIR is the near-infrared band, ρ Red It is a red light band, obtained through observations of optical remote sensing satellites such as MODIS, Landsat, and Sentinel-2. Satellite data with different temporal and spatial resolutions can be downloaded and used according to research needs. ini The choice of kernel function and manipulation of its parameters allow kNDVI ini Performs automatic and pixel-level adaptive stretching and ensures that all cases of relationship between the near-infrared and red bands are taken into account. This also makes kNDVI iniIt can cope with saturation effects, complex phenological cycles and seasonal changes, handle mixed pixel problems, and reduce uncertainty. Formula 23 is kNDVI ini The simplified calculation process has been proven to have achieved good application results in global research.

[0073] (2) Humidity index kSVMI

[0074] The kernel-driven moisture index kSVMI is a relatively simple remote sensing index based on the idea of ​​Gaussian kernel function, which can accurately invert regional soil moisture and be used in the process of carbon storage estimation. The kernel method (Kernel Function, KF) is used to calculate the square of the vector modulus and multiply it with the sign function of the corresponding SVMI (Index of Soil and Vegetation Moisture). This composite expression constitutes the new kernel soil moisture index kSVMI (Kernel-based optical remote sensing Index of Soil and Vegetation Moisture) formula proposed in this paper:

[0075] SVMI=3ρ Blue -ρ Red -ρ SWIR2 -ρ SWIR3 (Formula 25)

[0076]

[0077] In the formula, ρ Blue is the blue light band, ρ Red For the red light band, It is a short-wave infrared band, which is obtained through observations of optical remote sensing satellites such as MODIS, Landsat, Sentinel-2, etc. Satellite data with different temporal and spatial resolutions can be downloaded and used according to research needs. i can be 1, 2, and 3, representing short-wave infrared bands of different wavelength ranges.

[0078] (3) Carbon density index (CDI)

[0079] Carbon density index CDI is obtained by arithmetic average of kNDVI and kSVMI. iniRepresents the original kNDVI, whose value range is (-1, 1), kSVMI value range is (0, 1), and the theoretical value range of carbon density index should be (0, 1). kNDVI needs to be linearly stretched to make its value range (0, 1), and then arithmetic averaged with kSVMI. As a carbon density index CDI that couples the greenness index kNDVI and the humidity index kSVMI, its data type is a raster data type, which can indicate both greenness and humidity in surface quality.

[0080]

[0081] 3. Calculation of ecosystem carbon stocks taking into account land surface type and quality

[0082] We use the carbon density index to linearize the initial carbon storage to calculate the carbon storage that takes into account both type and quality, so that the carbon storage in the same land class reflects the quality difference. The carbon storage grid image of the kth part in a certain land class i is calculated using the improved method (IM) and To express, the specific calculation formula is shown in formula 30:

[0083]

[0084] Represents the maximum grid image of the theoretical carbon storage in the kth part of a certain land type i calculated by substituting the maximum annual average temperature and the maximum annual precipitation into the improved carbon density formula. Represents the minimum grid image of the theoretical carbon storage in the kth part of a certain land type i calculated by substituting the minimum annual average temperature and the minimum annual precipitation into the improved carbon density formula. CDI stands for carbon density index. For example, The maximum grid image of theoretical carbon storage in the aboveground part (k = above) of the forest land (i = 2) is obtained by substituting the highest annual average temperature and the highest annual precipitation into the improved carbon density formula. A carbon density index of 0 means that the pixel has the lowest quality among all the pixels in its land class, and its carbon storage is equal to the minimum carbon storage value under this land class; when the carbon density index is equal to 1, it means that the pixel has the highest quality among all the pixels in its land class, and the carbon storage takes the maximum carbon storage value corresponding to this land class. Similarly, carbon storage of different qualities under the same land class is linearized by the carbon density index within the numerical range of the maximum and minimum carbon storage corresponding to this land class. This method can more accurately indicate the difference in carbon storage in the study area to a certain extent without exceeding the theoretical maximum and minimum carbon storage calculated by the original method. After calculating the carbon storage grid image, the carbon storage of each part and the total carbon storage of each land class in the study area can be calculated by formula 31. The calculation principle is the same as formula 1.

[0085]

[0086] represents the total carbon storage value of a land type i in the study area calculated using the improved method, represents the aboveground carbon storage of land type i in the study area calculated using the improved method, and the same is true for the others. i Represents the area of ​​land type i.

[0087] Examples

[0088] Taking the Ningdong mining area as an example, the Sentiel-2 satellite observation data was used to improve the estimation of carbon reserves from 2017 to 2022. The results were evaluated from both qualitative and quantitative perspectives, and it was found that the improved results have been improved in both time correlation and spatial distribution accuracy.

[0089] The Ningdong Coal Base is geographically located at 106°22′31″-106°54′52″E, 37°4′41″-38°17′34″N. It is an approved key development zone and one of my country’s large-scale coal bases with an annual output of 100 million tons. It is also a demonstration zone for high-quality development, a high-tech industrial development zone and a chemical park. The mining area is dry and rainy, with rainfall mostly concentrated in July-September. The sunshine time is long and the evaporation is strong. The average annual precipitation is 255.2 mm. The vegetation types in the mining area are mainly sandy desert grasslands. The forest area is small, with many low mountains and hills. In addition to mining land, land for human activities such as residential buildings and cultivated land are relatively scarce. There are 8 small mining areas in the mining area, namely Hengcheng Mining Area, Hongdunzi Mining Area, Jijiajing Mining Area, Lingwu Mining Area, Majiatan Mining Area, Mengcheng Mining Area, Weizhou Mining Area and Yuanyang Lake Mining Area. The location diagram is shown below. Figure 2 shown.

[0090] Qualitative evaluation - spatial distribution comparison

[0091] The spatial distribution results Figure 3 and Figure 4 They represent the original carbon stock calculation results and the mining area carbon stock estimation results taking into account type and quality respectively.

[0092] In terms of spatial distribution, the improved carbon stock calculation results are more detailed in space. The carbon stocks in the same land type with different surface qualities show differences, realizing the carbon stock calculation that takes into account both surface type and quality.

[0093] Quantitative evaluation - comparison of aboveground carbon storage

[0094] (1) Inter-class correlation comparison

[0095] Leaf Area Index (LAI) remote sensing products were selected for the evaluation of aboveground carbon storage. LAI is a measure of the ratio of plant leaf area to land surface area. It is usually used to evaluate the growth status and biological productivity of vegetation. It is one of the important parameters for global carbon cycle, climate change research and agricultural production. A higher LAI means that plants have more leaves for photosynthesis. Plants absorb carbon dioxide from the atmosphere through photosynthesis and convert it into organic matter, thereby increasing the carbon input of the ecosystem, which is conducive to the increase of carbon storage. In the experiment, MODIS / 061 / MOD15A2H and MODIS / 061 / MYD15A2H were selected as LAI products, with a spatial resolution of 500m and a composite data set of 8 days. By calculating the mean carbon storage obtained by the original method, the mean carbon storage obtained by the improved method and the mean LAI product in each land category from 2017 to 2022, the correlation between carbon storage and products between land categories before (OM) and after (IM) was obtained to judge the reliability of the improved method. The correlation is shown in Tables 2 and 3 below:

[0096] Table 2 Comparison of correlation between aboveground carbon storage and MOD15A2H product types

[0097]

[0098] Table 3 Comparison of correlation between aboveground carbon storage and MYD15A2H product categories

[0099]

[0100] The inter-class correlations between carbon storage and LAI products before and after the improvement were high, and the Pearson correlation coefficients mostly remained near 0.4. This result firstly showed that the carbon storage calculated based on the InVEST carbon storage and carbon sequestration model had good consistency with the LAI products, which could reflect the regional carbon storage level to a certain extent. Secondly, the inter-class correlations between carbon storage and LAI obtained by the improved method were improved during the study year (green represents better indicators, red represents worse indicators, and the calculation results of each indicator are retained to 3 decimal places), and the improvement was between 1.50% and 5.00%. The verification results of the above inter-class correlations further demonstrated that the carbon storage obtained based on the InVEST model and improved on its basis is effective and performs better.

[0101] (2) Intra-class correlation comparison

[0102] The verification of intra-class correlation is reflected by calculating the pixel-by-pixel temporal correlation between carbon storage and LAI products from 2017 to 2022. Since the carbon storage calculated by the method before the improvement does not differ within the class, the carbon storage values ​​in the same land class are the same. Therefore, when calculating the temporal correlation, the standard deviation of carbon storage is 0, and Pearson correlation does not exist. This is also the defect of the original method. Taking the MYD15A2H product as an example, the temporal correlation between the improved carbon storage and the product can be obtained, and its spatial distribution map is shown in the figure below. Figure 5 As shown, the intra-class correlation analysis statistics are shown in Table 4 below.

[0103] Table 4 Statistical table of correlation between aboveground carbon storage and MYD15A2H product category

[0104]

[0105] Qualitative results Figure 5 It shows that in most areas of the Ningdong mining area, the time correlation between the improved carbon reserves and each product is high and positively correlated; the time correlation in the Yuanyang Lake mining area and its surrounding areas, that is, the northern part of the Ningdong mining area, is partially low, but the range is not large.

[0106] The intra-class correlation results in Table 4 show that the forest land has a low correlation and the grassland has a high correlation. The forest land accounts for a very small proportion in the Ningdong mining area, less than 2%, and the number of samples is small. The larger the LAI, the more leaves it has, which usually means the older the tree is. The growth of old trees is slower than that of new trees, and the less organic matter they provide, the weaker their carbon sequestration capacity. Most of the trees in the forest land are perennial trees, that is, there are more old trees, so the carbon storage in the forest land is not positively correlated with the LAI. This is another important reason why the temporal correlation between the carbon storage in the forest land and the LAI product is negative. As the land type with the largest area, the high temporal correlation of the grassland part further illustrates that the carbon storage calculated using the improved method has a high consistency with the LAI product in time, and has the ability and feasibility to continuously draw carbon storage in time.

[0107] Quantitative evaluation - comparison of carbon storage in other parts

[0108] The other carbon stocks are closely related to underground and soil carbon, so the Soil organic carbon stock in the SoilGrids 2.0 product is selected for comparison and verification with the calculated carbon stocks of other parts. SoilGrids 2.0 is a global digital soil map system that uses more than 230,000 soil profile observation data and a series of environmental covariates (including environmental information such as climate, land cover and topography) in the WoSIS database, combined with machine learning methods to map the spatial distribution of global soil physical and chemical properties, with a spatial resolution of 250m. The above are the remote sensing products used in the current carbon stock calculation and its improvement process. After collection and download, they were pre-processed by reprojection, resampling and clipping, and then added to the comparative verification experiment. The Soil organic carbon stock is regarded as the true value of underground or soil carbon storage, and the results before and after improvement are evaluated using indicators including RMSE, MAE, bias, Pearson correlation coefficient and WIA. The accuracy evaluation indicators involved are formula 32-formula 36.

[0109]

[0110] where x i ,y i represent the observed values ​​of two variables, and The formulas for the mean values ​​of the representative samples are and In this study, the Pearson correlation coefficient is used to verify the trend of the calculation results and the remote sensing products, with a range of [-1, 1].

[0111]

[0112] P i Represents the predicted value of the i-th observation, that is, the calculated value, O i represents the i-th observation, that is, the product value or the value considered to be the true value, and n represents the number of samples.

[0113]

[0114] Where n is the number of samples, y i is the true value of the i-th sample, is the predicted value of the i-th sample. MAE calculates the absolute difference between the predicted value and the actual value, and averages these differences to obtain the overall error size.

[0115]

[0116] The WIA index is a standardized metric index for the degree of model prediction error proposed by Wilmot. Among them, O represents the observed value, P represents the predicted value, Represents the average observation value. WIA takes values ​​between 0 and 1. Represents the ratio of the mean square error to the potential error. 1 means that the estimated value matches the actual value perfectly, and 0 means that the estimated value does not match the actual value at all. The calculation method of the average observation value in the above formula is:

[0117]

[0118] Table 5 Comparison of different combinations of other parts with SoilGrids 2.0 products

[0119]

[0120] In order to explore the consistency between the different combination results of the carbon stock of other parts and the soil product, three combinations were set. In addition to the soil part, the Soil organic carbon stock also includes the organic carbon part. The verification results in Table 5 show that whether it is the sum of the soil carbon stock C_soil and the underground carbon stock C_below, or the sum of the soil carbon stock C_soil and the dead organic matter carbon stock C_dead, or the sum of the three parts, the evaluation indicators all show that the carbon stock (IM) calculated using the improved method is better (green represents better indicators, red represents worse indicators, and the calculation results of each indicator are retained to 3 decimal places), indicating that the method of calculating the carbon stock of the mining area by taking into account the type and quality has a certain applicability in the underground and soil parts, and can obtain more accurate results than the original method on the basis of fully expressing spatial heterogeneity. The above example verification experiment proves the feasibility and advantages of the improved method from both qualitative and quantitative perspectives.

[0121] Unless otherwise defined, all technical and / or scientific terms used in this application have the same meaning as commonly understood by one of ordinary skill in the art to which the invention relates. The materials, methods and examples mentioned in this application are illustrative only and not limiting.

[0122] Although the present invention has been described in conjunction with specific embodiments, within the spirit of the invention of this application, those skilled in the art may make appropriate substitutions, modifications and changes, and such substitutions, modifications and changes still fall within the scope of protection of this application.

Claims

1. A remote sensing calculation method for ecosystem carbon storage that takes into account both surface type and quality, comprising: The average carbon density table, maximum carbon density table, minimum carbon density table and land use data of the study area were imported into the carbon storage and sequestration module of InVEST software, and the theoretical original carbon storage raster image of each type of land in the study area was calculated through its built-in algorithm. Theoretical maximum original carbon storage raster image and theoretical minimum original carbon storage raster images The resolution of the imagery is consistent with the land use data; use To calculate the carbon storage grid image of the kth part in a certain land type i Among them, k represents a part of carbon storage, k represents the total carbon storage when it is total, k represents the aboveground carbon storage when it is above, k represents the underground carbon storage when it is below, k represents the soil carbon storage when it is soil, and k represents the dead organic carbon storage when it is dead; i represents the number of each type of land, 1 for cultivated land, 2 for forest land, 3 for grassland, 4 for water area, 5 for building land, and 6 for unused land; CDI is the carbon density index, which is used to indicate the surface quality; is the aboveground carbon stock of land type i calculated using formula 30; is the underground carbon stock of land type i calculated using formula 30; is the carbon stock in the soil of land type i calculated using formula 30; is the dead organic carbon stock of land type i calculated using formula 30; Then use To calculate the total carbon storage of land type i in the study area, A i Represents the area of ​​land type i; the sum of the total carbon storage of each land type is the total carbon storage of the study area.

2. The method according to claim 1, characterized in that: CDI according to Calculation is performed; kNDVI is the greenness index; kSVMI is the humidity index; in, Among them, tanh is the hyperbolic tangent function, sgn is the sign function, and ρ NIR is the near-infrared band, ρ Red It is the red light band; In the formula, ρ Blue is the blue light band, ρ Red For the red light band, is the short-wave infrared band, and j can be 1, 2 and 3, representing short-wave infrared bands of different wavelength ranges.

3. The method according to claim 1, characterized in that: in, is the soil carbon density obtained based on the annual maximum precipitation HAP, and are the biomass carbon density obtained according to the annual maximum temperature HAT and the annual maximum precipitation HAP; C′ i and C″ i Data represent the study area and the whole country, respectively; is the maximum biomass carbon density correction factor; is the maximum soil carbon density correction factor Use the maximum biomass carbon density correction factor Multiply them with the aboveground carbon density, underground carbon density and dead organic matter carbon density columns in the national carbon density data, and use the maximum soil carbon density correction factor Multiply it with the soil carbon density column in the national carbon density data to obtain the maximum carbon density table of the study area.

4. The method according to claim 1, characterized in that: in, is the soil carbon density obtained based on the annual minimum precipitation LAP, and are the biomass carbon density obtained according to the annual minimum temperature LAT and the annual minimum precipitation LAP; C′ i and C″ i Data represent the study area and the whole country, respectively; is the minimum biomass carbon density correction factor, is the minimum soil carbon density correction factor; Use minimum biomass carbon density correction factor Multiply them with the aboveground carbon density, underground carbon density and dead organic matter carbon density columns in the national carbon density data, and use the minimum soil carbon density correction factor Multiplying it with the soil carbon density column in the national carbon density data, the minimum carbon density table of the study area was obtained.

5. The method according to claim 1, characterized in that: C SP =3.3968×MAP+3996.1 (R 2 =0.11) (Formula 2); C BP =6.798×e 0.0054×MAP (R 2 =0.70) (Formula 3); C BT =28×MAT+398 (R 2 =0.47, P<0.01) (Formula 4); K B =K BT ×K BP (Formula 7); Among them, C SP is the soil carbon density obtained based on the average annual precipitation MAP, C BT and C BP are the biomass carbon density obtained based on the average annual temperature MAT and the average annual precipitation MAP; C′ i and C″ i represent the data of the study area and the whole country respectively. When the biomass carbon density is improved according to precipitation, i = BP; when the biomass carbon density is improved according to temperature, i = BT; when the soil carbon density is improved according to precipitation, i = SP; K B is the biomass carbon density correction factor; K S is the soil carbon density correction factor; the biomass carbon density correction factor K is used B Multiply them with the aboveground carbon density, underground carbon density and dead organic matter carbon density columns in the national carbon density data, and use the soil carbon density correction factor K S Multiplying it with the soil carbon density column in the national carbon density data, we get the average carbon density table for the study area.

Images