Soil carbon loss remote sensing image monitoring method and system in hydraulic erosion area

Abstract

The application discloses a method and system for monitoring soil carbon loss in a water erosion area through remote sensing images, and belongs to the field of remote sensing image processing.The method comprises the following steps: extracting spectral characteristics by using a hyperspectral image, adopting a PLSR and SVM algorithm to construct a model, and inverting soil organic carbon content; superimposing the carbon content and an erosion intensity classification map, calculating a carbon loss rate, and establishing a coupling relationship; combining terrain parameters to construct a source-sink index, dividing an erosion carbon output area, a migration transition area, and a deposited carbon enrichment area; and comparing carbon content changes during a snow melting period and before and after a rainy season, and accounting for seasonal carbon loss.

Description

Technical Field

[0001] This invention relates to the fields of remote sensing image processing and soil and water conservation monitoring technology, specifically to a method and system for monitoring soil carbon loss in water erosion zones using remote sensing images. Background Technology

[0002] Water erosion is a significant driver of soil organic carbon loss, particularly in semi-arid and high-altitude mountainous areas. Due to the combined effects of seasonal freeze-thaw cycles and concentrated rainfall, soil erosion is accompanied by substantial lateral migration of organic carbon. Statistics show that billions of tons of soil organic carbon are lost globally each year due to water erosion, a considerable proportion of which occurs in ecologically fragile areas such as high-altitude mountainous regions. As the largest carbon sink component of terrestrial ecosystems, soil organic carbon loss not only leads to decreased soil fertility and ecological degradation but may also exacerbate the greenhouse effect by releasing carbon dioxide through mineralization. Therefore, accurately monitoring the spatiotemporal dynamics of soil carbon loss in water erosion areas is of significant scientific and practical value for regional carbon cycle research and assessment of carbon sequestration benefits in soil and water conservation.

[0003] The migration and transformation of soil organic carbon during erosion is a complex physicochemical process involving multiple stages such as carbon stripping, transport, deposition, and mineralization. In erosion-affected areas, topsoil rich in organic carbon is stripped by runoff and transported downhill, leading to a decrease in in-situ carbon content. In depositional areas, sediment carrying organic carbon is deposited and enriched, resulting in a relative increase in carbon content in the region. This erosion-driven lateral redistribution of carbon profoundly affects the carbon balance at the watershed scale. Especially in high-altitude and cold mountainous areas, the spring snowmelt season and the summer rainy season are the two main active periods for hydraulic erosion, with snowmelt and rainwater runoff driving carbon loss processes of varying intensities and spatial patterns.

[0004] Remote sensing technology, due to its advantages of acquiring surface information over a large area, rapidly, and non-destructively, has become an important means of monitoring soil erosion and carbon cycle. Existing remote sensing monitoring technologies for soil erosion mainly focus on the extraction and analysis of erosion morphology characteristics. For example, Chinese patent application CN116934753A discloses a soil and water conservation monitoring method based on remote sensing images. This method extracts the edge lines of erosion gullies by edge detection on high-resolution remote sensing images, constructs the gully span ratio based on the direction and width changes of the gullies, and constructs a unilateral gully subsidence index by combining the correspondence between gully depth and image grayscale values. Finally, it calculates the soil erosion significance to assess the degree of soil erosion. This method improves the image segmentation accuracy of soil erosion areas to a certain extent, providing technical support for the evaluation of soil and water conservation effects.

[0005] In the field of soil organic carbon remote sensing inversion, the development of hyperspectral remote sensing technology has provided a new approach for the areal estimation of soil carbon content. Soil organic matter exhibits characteristic absorption properties in the visible and near-infrared bands. Utilizing this spectrally sensitive response characteristic, rapid estimation of soil carbon content can be achieved by establishing a spectral inversion model. Machine learning algorithms such as partial least squares regression and support vector machines have been widely applied to spectral modeling of soil carbon content, achieving good predictive results. However, existing research mostly focuses on soil carbon content inversion in farmland or plains areas, with relatively few studies on water erosion areas, especially high-altitude and cold mountainous regions, and few systematic methods combining carbon content inversion with erosion process analysis.

[0006] However, the existing technologies still have the following shortcomings. First, existing methods mainly focus on the physical morphological characteristics of soil erosion, such as the geometric parameters of erosion gullies, but fail to establish a quantitative correlation between the erosion process and soil organic carbon loss, making it difficult to directly obtain spatial distribution information on carbon loss. Second, existing technologies lack a systematic analysis of the coupling relationship between erosion intensity and carbon loss, failing to reveal the differences and changing patterns of carbon content in areas with different erosion levels. This lack of coupling makes it difficult to understand the degree of erosion's impact on the carbon pool from a mechanistic perspective. Third, existing methods fail to identify the source-sink spatial pattern of organic carbon in the study area, making it difficult to track the spatial redistribution process of organic carbon within the watershed. This source-sink pattern information is of great guiding significance for formulating differentiated soil and water conservation measures. Fourth, existing technologies lack the ability to conduct temporal comparative analysis of snowmelt and rainy seasons in high-altitude and cold mountainous areas, failing to distinguish the contribution of different seasonal erosion processes to carbon loss, which limits a deeper understanding of the spatiotemporal dynamics of carbon loss.

[0007] Therefore, there is an urgent need to develop a remote sensing monitoring method and system for soil carbon loss in water erosion areas that can comprehensively utilize hyperspectral remote sensing data to invert soil carbon content, establish an erosion-carbon loss coupling analysis mechanism, identify carbon loss source and sink patterns, and achieve seasonal time-series monitoring, so as to fill the gap in the existing technology for quantitative monitoring of carbon loss in erosion areas. Summary of the Invention

[0008] To address the shortcomings of existing technologies, this invention provides a remote sensing image monitoring method and system for soil carbon loss in water erosion areas, aiming to achieve quantitative monitoring and spatial analysis of soil organic carbon loss in water erosion areas.

[0009] This invention provides a method for monitoring soil carbon loss in water erosion zones using remote sensing images, comprising the following steps:

[0010] Step S1, Spectral Inversion Step: Acquire hyperspectral remote sensing images of the study area, preprocess the hyperspectral remote sensing images to obtain spectral reflectance data, extract spectral features in the visible-near-infrared bands, identify bare soil pixels, combine with measured soil organic carbon content data from ground sampling points, construct a spectral inversion model using partial least squares regression and support vector machine algorithms, perform areal estimation of soil organic carbon content in the surface layer of bare soil, and obtain soil carbon content inversion results;

[0011] Step S2, Coupling Analysis Step: Obtain the erosion intensity classification map of the study area, perform spatial overlay analysis between the soil carbon content inversion results and the erosion intensity classification map, calculate the statistical characteristic value of carbon content and the carbon content difference value of different erosion level areas, calculate the carbon loss rate based on the erosion intensity and carbon content difference value, and establish the erosion-carbon loss coupling relationship.

[0012] Step S3, source-sink identification step: Obtain digital elevation model data of the study area, extract topographic parameters such as topographic slope, topographic humidity index and runoff accumulation, construct source-sink pattern discrimination index in combination with erosion and deposition characteristics, divide the study area into erosion carbon output area, migration transition area and deposition carbon enrichment area, and identify the spatial redistribution pattern of organic carbon in the watershed.

[0013] Step S4, temporal monitoring steps: acquire multi-temporal hyperspectral remote sensing images before and after the snowmelt period and before and after the rainy season, select images of four key time phases, execute steps S1 to S3 respectively, compare the spatial distribution changes of soil carbon content in different time phases, calculate the carbon loss due to snowmelt erosion and the carbon loss due to rainfall erosion respectively, and output the spatiotemporal dynamic monitoring results of carbon loss.

[0014] Another aspect of the present invention provides a remote sensing image monitoring system for soil carbon loss in water erosion zones, comprising:

[0015] The data acquisition module is used to acquire hyperspectral remote sensing images, erosion intensity classification maps, digital elevation model data, and measured data of soil organic carbon content at ground sampling points in the study area.

[0016] The spectral inversion module is used to preprocess hyperspectral remote sensing images to obtain spectral reflectance data, extract spectral features in the visible-near-infrared bands, construct a spectral inversion model using partial least squares regression and support vector machine algorithms, and obtain soil carbon content inversion results.

[0017] The coupling analysis module is used to perform spatial overlay analysis of soil carbon content inversion results and erosion intensity classification map, calculate the statistical characteristic value of carbon content and the carbon content difference value of different erosion level areas, and calculate the carbon loss rate based on erosion intensity and carbon content difference value.

[0018] The source-sink identification module is used to extract topographic parameters and combine them with erosion and sedimentation characteristics to construct a source-sink pattern discrimination index, which divides the study area into an erosion carbon output zone, a migration transition zone, and a sedimentary carbon enrichment zone.

[0019] The time-series monitoring module is used to acquire multi-temporal hyperspectral remote sensing images, compare the spatial distribution changes of soil carbon content in different time phases, and calculate the carbon loss due to snowmelt erosion and rainfall erosion, respectively.

[0020] The beneficial effects of this invention are as follows: This invention achieves areal estimation of soil organic carbon content through spectral inversion technology of hyperspectral remote sensing images, overcoming the spatial coverage limitations of traditional point sampling; by constructing an erosion-carbon loss coupling analysis module, a quantitative relationship between erosion intensity and carbon loss is established, revealing the differences in carbon content in areas with different erosion levels; by designing a source-sink pattern identification algorithm, the spatial redistribution process of organic carbon within the watershed is tracked; and by using a time-series monitoring module, the contributions of snowmelt erosion and rainfall erosion to carbon loss are distinguished, providing an effective remote sensing monitoring method for carbon cycle research in high-altitude and cold mountainous areas and for assessing the carbon sink benefits of soil and water conservation. Attached Figure Description

[0021] Figure 1 This is a flowchart of a remote sensing image monitoring method for soil carbon loss in water erosion areas provided in an embodiment of the present invention.

[0022] Figure 2 This is an architecture diagram of the remote sensing image monitoring system for soil carbon loss in water erosion areas provided in an embodiment of the present invention. Detailed Implementation

[0023] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only a part of the embodiments of the present invention, and not all of them. All other embodiments obtained by those skilled in the art based on the embodiments of the present invention without creative effort are within the scope of protection of the present invention.

[0024] See Figure 1 This invention provides a remote sensing image monitoring method for soil carbon loss in water erosion areas. This method is applicable to dynamic monitoring of soil carbon loss in water erosion-sensitive areas such as semi-arid and high-altitude mountainous regions. Its core idea lies in organically combining hyperspectral remote sensing technology, erosion process analysis, and carbon cycle research to form a four-in-one monitoring technology system encompassing spectral inversion, coupled analysis, source-sink identification, and time-series monitoring. The overall method includes the following steps:

[0025] Step S1, spectral inversion step.

[0026] The main purpose of this step is to utilize the detailed spectral information from hyperspectral remote sensing images to invert the spatial distribution of soil organic carbon content. In one embodiment of the present invention, this step specifically includes the following sub-steps.

[0027] First, acquire hyperspectral remote sensing images of the study area. Preferably, GF-5 hyperspectral satellite imagery, HJ-1A / 1B hyperspectral imagery, or image data acquired by an airborne hyperspectral imaging system can be used. The spectral range of the hyperspectral remote sensing images should cover 400nm to 2500nm, with a spectral resolution of not less than 10nm. The spatial resolution is selected according to the scale of the study area; for watershed-scale studies, a resolution of 10m to 30m is recommended. When acquiring images, clear, cloudless weather conditions should be selected. Images taken during the period of surface exposure (such as after snowmelt in spring or after autumn harvest) are more conducive to the extraction of soil spectral information. In one embodiment of the present invention, a typical alpine grassland watershed in the eastern section of the Qilian Mountains in Qinghai Province is used as the study area. The elevation range of this watershed is 2800m to 4200m, with an average annual precipitation of about 450mm, of which more than 60% is concentrated between July and September. The average annual temperature is about -1.5 degrees Celsius, and permafrost is widely distributed. The study area is dominated by snowmelt erosion and rainfall erosion, with erosion intensity ranging from slight to severe, making it an ideal area for conducting carbon loss monitoring research in erosion zones.

[0028] Secondly, hyperspectral remote sensing images are preprocessed to obtain spectral reflectance data. The preprocessing process includes three main steps: radiometric correction, atmospheric correction, and geometric correction. Radiometric correction converts the digital quantization values ​​of the image into radiance values. Its basic principle is to use the sensor's calibration coefficients to perform a linear transformation on the original DN values. The specific formula is as follows: ,in Radiance, expressed in W / (m²·sr·μm). and These are the gain coefficient and offset, which can be obtained from the image header file or calibration parameter file. Atmospheric correction uses the FLAASH model or 6S model to eliminate the effects of atmospheric scattering and absorption, converting radiance values ​​into surface reflectance. In one embodiment of the invention, the FLAASH model is preferably used for atmospheric correction. This model is based on the MODTRAN radiative transfer code and can handle atmospheric conditions in high-altitude areas well. The correction process requires inputting parameters such as image acquisition time, sensor type, atmospheric model (usually subarctic winter or subarctic summer for high-altitude mountainous areas), and aerosol model (usually rural type). Geometric correction uses a combination of ground control points and digital elevation models to ensure accurate registration of the image with the geographic coordinate system, controlling the positional error to within 0.5 pixels. For mountainous images, orthorectification is also required to eliminate geometric distortion caused by terrain undulations. Orthorectification is achieved using a rational function model or a rigorous physical model combined with DEM data.

[0029] Furthermore, spectral characteristics in the visible-near-infrared band are extracted. Soil organic matter exhibits characteristic absorption properties in the visible-near-infrared band, mainly manifested in the following sensitive band ranges: the 480nm to 520nm band reflects the absorption characteristics of soil organic matter to blue light; the higher the organic matter content, the lower the reflectivity in this band; the 580nm to 680nm band reflects the complexation effect between iron oxides and organic matter; organic matter affects the reflectivity in this band by complexing with iron oxides; absorption peaks near 1400nm, 1900nm, and 2200nm are associated with soil moisture and clay minerals, indirectly reflecting the occurrence environment of organic carbon. Among them, the absorption near 1400nm and 1900nm is mainly caused by the OH bond vibration of water molecules, and the absorption near 2200nm is related to the Al-OH bond vibration in clay minerals. In one embodiment of the present invention, the original spectrum is enhanced using the continuum removal method and first-order differential transformation to highlight the characteristic absorption information of organic carbon. The continuum removal method enhances the contrast of absorption features by constructing the convex envelope of the spectral curve and calculating the normalized reflectance. The calculation formula is as follows: ,in The reflectance after removing the continuum is . Original reflectivity This represents the value at the wavelength corresponding to the envelope; the first-order differential transform eliminates the influence of baseline shift by calculating the difference in reflectance between adjacent bands, making the spectral characteristics more obvious. Its calculation formula is: ,in wavelength The first derivative value at the specified point. After spectral enhancement processing, the characteristic bands with the highest correlation to soil organic carbon content are selected as model input variables. In one embodiment of the present invention, the Pearson correlation coefficient between each band and the measured carbon content is calculated, and bands with an absolute value of correlation coefficient greater than 0.4 and passing the significance test are selected as candidate characteristic bands. Then, the optimal combination of characteristic bands is further selected using a continuous projection algorithm or a competitive adaptive reweighted sampling algorithm.

[0030] Simultaneously, measured data on soil organic carbon content from ground sampling points are needed for model establishment and validation. The layout of sampling points should consider both spatial representativeness and gradient coverage of erosion intensity. Preferably, the sampling density is no less than 3 sampling points per square kilometer, with a sampling depth of 0 to 20 cm in the topsoil. Soil organic carbon content is determined using the potassium dichromate oxidation-external heating method or elemental analysis, with results recorded in g / kg. During sampling, auxiliary information such as GPS coordinates, topographic location, and vegetation cover should be recorded simultaneously to facilitate subsequent spatial analysis and model correction.

[0031] Then, a spectral inversion model is constructed using partial least squares regression and support vector machine algorithms. In one embodiment of this invention, the model construction employs a strategy of parallel modeling of two algorithms and selection of the optimal one. The partial least squares regression (PLSR) model is suitable for cases where there is multicollinearity among spectral variables. Its basic principle is to reduce the dimensionality of independent variables by extracting principal components while retaining the information most strongly correlated with the dependent variable. The support vector machine (SVM) model maps the low-dimensional feature space to a high-dimensional feature space through a kernel function, and finds the optimal classification hyperplane in the high-dimensional space, making it suitable for handling nonlinear relationships.

[0032] For the partial least squares regression model, the inversion formula for soil organic carbon content can be expressed as:

[0033] ,

[0034] in: The value is the soil organic carbon content inversion value, with the unit being g / kg. Its value ranges from 1 to 80 g / kg, reflecting the actual variation range of soil organic carbon content in the study area. The intercept term represents the baseline carbon content value when all spectral variables are zero. Its value is obtained by fitting the modeling samples using the least squares method. For the first The regression coefficients of the spectral characteristic variables are dimensionless, and their positive and negative signs and absolute values ​​reflect the direction and intensity of the correlation between the band and the organic carbon content, respectively. For the first The values ​​of the spectral characteristic variables can be the original reflectance, the continuum removed value, or the first derivative value, and the values ​​range from 0 to 1. The number of feature bands used in modeling is typically between 5 and 15, with the specific number determined through cross-validation. The rationale for selecting 5 to 15 feature bands is that too few bands will result in the loss of effective spectral information, while too many bands will introduce noise and lead to overfitting. The optimal number of bands is determined through 10-fold cross-validation, ensuring the coefficient of determination of the validation set is optimal. maximize.

[0035] For the support vector machine model, one embodiment of the present invention employs a radial basis function (RBF), with kernel function parameters... and penalty parameters The optimal value is determined by combining grid search with cross-validation. Preferably, The search range is from 0.001 to 1000. The search range is from 0.1 to 10000, all divided using logarithmic intervals. The input to the support vector machine model is the standardized spectral values ​​of the feature bands, and the output is the predicted soil organic carbon content.

[0036] Model validation employs leave-one-out cross-validation or independent validation set validation, with evaluation metrics including the coefficient of determination. Root mean square error (RMSE) and relative analysis error (RPD). In one embodiment of the present invention, when A value greater than 0.65, an RMSE less than 3.5 g / kg, and an RPD greater than 1.8 indicate that the model has good predictive ability and can be used for spatial inversion. If the validation accuracy of the two models is not significantly different, the partial least squares regression model is preferred because it has clearer physical meaning and better stability; if the accuracy of the support vector machine model is significantly higher than that of the partial least squares regression model, the support vector machine model is used for inversion.

[0037] Finally, the constructed spectral inversion model was applied to the hyperspectral imagery of the entire region to perform areal estimation of soil organic carbon content in the exposed soil surface layer, obtaining the soil carbon content inversion results. During the inversion process, soil pixel identification is required first, removing non-bare soil pixels such as those in vegetation-covered areas, water bodies, and built-up areas. In one embodiment of this invention, a threshold of Normalized Difference Vegetation Index (NDVI) less than 0.2 is used to identify bare soil pixels. This threshold is selected based on the fact that when NDVI is less than 0.2, vegetation cover is usually below 15%, soil spectral information dominates, and interference with carbon content inversion is minimal. The inversion results are stored in the form of a raster layer, with each pixel recording the estimated soil organic carbon content at its corresponding location, forming a spatial distribution map of soil carbon content in the study area.

[0038] Step S2, coupling analysis step.

[0039] The main objective of this step is to establish a quantitative correlation between soil erosion intensity and organic carbon loss, revealing the driving mechanism of carbon loss by the erosion process. In one embodiment of the present invention, this step specifically includes the following sub-steps.

[0040] First, an erosion intensity classification map of the study area is obtained. This map can be calculated using the Modified Universal Soil Loss Equation (RUSLE) model or existing soil and water conservation survey data. In one embodiment of this invention, the calculation formula for the RUSLE model is as follows: ,in This represents the average annual soil erosion. As the erosivity factor of rainfall, As a soil erodibility factor, For terrain factors, As a factor of vegetation cover and management, This is a factor for soil and water conservation measures. According to the "Classification and Grading Standard for Soil Erosion" (SL190-2007), the calculated soil erosion modulus is divided into six levels: slight erosion (less than 500 t / (km²·a)), mild erosion (500 to 2500 t / (km²·a)), moderate erosion (2500 to 5000 t / (km²·a)), severe erosion (5000 to 8000 t / (km²·a)), extremely severe erosion (8000 to 15000 t / (km²·a)), and extreme erosion (greater than 15000 t / (km²·a)).

[0041] Secondly, the soil carbon content inversion results were spatially overlaid with the erosion intensity classification map. The spatial overlay operation was performed in the Geographic Information System (GIS) platform, matching the soil carbon content raster layer with the erosion intensity classification layer at the pixel level, ensuring that each raster unit simultaneously possesses both carbon content and erosion level attributes. Before overlaying, it was ensured that the two layers had the same spatial reference system and pixel size; if the resolutions were inconsistent, they were uniformly resampled to a coarser resolution.

[0042] Then, the statistical characteristic values ​​and carbon content difference values ​​of different erosion levels are calculated. For each erosion level region, the mean, standard deviation, maximum, minimum, and coefficient of variation of carbon content for all pixels within the region are statistically analyzed. In one embodiment of the present invention, the calculation of carbon content difference values ​​is based on a slightly eroded area, and the formula for calculating the carbon content difference values ​​between each erosion level region and the slightly eroded area is as follows:

[0043] ,

[0044] in: For the first The carbon content difference value of the erosion zone is expressed in g / kg. A positive value indicates that the carbon content of the erosion zone is lower than that of the benchmark zone, and a negative value indicates that the carbon content is higher than that of the benchmark zone. The average soil organic carbon content in the baseline area (slightly eroded area) is expressed in g / kg. For the first Mean soil organic carbon content in graded erosion zones, in g / kg; subscript The values ​​range from 1 to 6, corresponding to six erosion levels from slight to severe. The carbon content difference reflects the cumulative effect of carbon loss caused by the erosion process; the larger the difference value, the stronger the reduction effect of erosion on the carbon pool.

[0045] Furthermore, the carbon loss rate is calculated based on the difference between erosion intensity and carbon content. This invention proposes the concept of an erosion-carbon loss coupling coefficient to characterize the correlation strength between erosion intensity and carbon loss. The formula for calculating the coupling coefficient is:

[0046] ,

[0047] in: The erosion-carbon loss coupling coefficient ranges from 0 to 1 and is dimensionless. The closer the value is to 1, the stronger the coupling relationship between erosion and carbon loss. For the first The area weight of the graded erosion zone is equal to the ratio of the area of ​​that grade to the total area of ​​the study area, and is dimensionless. For the first The carbon content difference between graded erosion zones, in g / kg; The maximum value of the carbon content difference across all erosion grades, expressed in g / kg, is used for normalization. For the first The annual average soil erosion modulus of the grade is expressed in t / (ha·a), and the median value of the range of erosion modulus for that grade is taken. This represents the upper limit of erosion intensity for severe erosion levels, expressed in t / (ha·a), and is used for normalization. This coupling coefficient comprehensively considers the spatial distribution characteristics of erosion intensity and the response of carbon content, avoiding the one-sidedness that may result from a single indicator.

[0048] Based on the above analysis, the formula for estimating the carbon loss rate is:

[0049] ,

[0050] in: The carbon loss rate, expressed in g / (m²·a), represents the net loss of organic carbon per unit area per year due to erosion. is the erosion-carbon loss coupling coefficient, which is dimensionless and reflects the effective contribution rate of erosion to carbon loss. For the first The annual average soil erosion modulus of the grade, expressed in t / (ha·a); This refers to soil bulk density, expressed in g / cm³, with a typical range of 1.1 to 1.6 g / cm³. Specific values ​​can be obtained from actual measurements or by consulting regional soil databases. The effective erosion depth is expressed in cm and represents the thickness of the soil layer that is mainly affected by the erosion process. For surface erosion and gully erosion, it is usually taken as 5 to 10 cm, and for gully erosion, it can be taken as 10 to 30 cm. For the first The average organic carbon content of soil in graded erosion zones is expressed in g / kg; a constant of 1000 is used for unit conversion. The physical meaning of this formula is as follows: the erosion modulus reflects the amount of soil loss per unit area; soil bulk density converts volume into mass; effective depth defines the range of soil layers involved in erosion; carbon content represents the concentration of organic carbon in the eroded soil; and the coupling coefficient corrects for the actual proportion of carbon loss in the eroded soil.

[0051] Finally, an erosion-carbon loss coupling model is established, outputting a spatial distribution map of carbon loss rates for each erosion level and a coupling relationship analysis report. In one embodiment of the invention, the coupling relationship can be fitted using linear regression, exponential functions, or power functions, and the function form with the highest goodness of fit is selected as the final coupling relationship expression.

[0052] Step S3, source and sink identification step.

[0053] The main purpose of this step is to identify the spatial redistribution pattern of organic carbon within the study area and determine the distribution of carbon loss source and sink regions. In one embodiment of the present invention, this step specifically includes the following sub-steps.

[0054] First, digital elevation model (DEM) data for the study area is acquired, and topographic parameters are extracted. The spatial resolution of the DEM should match or be finer than that of the hyperspectral imagery, preferably using DEM data with a resolution of 10m to 30m. Based on the DEM, the following topographic parameters are extracted: terrain slope. The topographic moisture index (TWI) indicates the degree of surface inclination, measured in degrees, and ranges from 0 to 90 degrees. It is dimensionless and reflects the potential for surface moisture accumulation; the calculation formula is... ,in The catchment area is the area per unit contour line length. The slope is the tangent of the slope; the runoff accumulation reflects the degree of runoff accumulation at the current location from upstream; the slope index reflects the relative height of the current location on the slope.

[0055] Secondly, a source-sink pattern discrimination index is constructed by combining erosion and deposition characteristics. This invention proposes a source-sink index (SI) that comprehensively considers topographic location and erosion and deposition characteristics, and its calculation formula is as follows:

[0056] ,

[0057] in: The source-sink index has a value range of -1 to 1. A positive value indicates a carbon exporting region (source), a negative value indicates a carbon enrichment region (sink), and a value close to 0 indicates a transitional region. The terrain slope of the current pixel, in degrees; The critical slope threshold is the minimum slope at which erosion occurs. For grasslands in high-altitude and cold mountainous areas, it is usually taken as 5 to 8 degrees. The basis for determining this threshold is that when the slope is less than the critical value, the surface runoff velocity is insufficient to move soil particles, and erosion is unlikely to occur. The maximum slope value in the study area is expressed in degrees. The topographic humidity index for the current pixel is dimensionless. and These are the maximum and minimum values ​​of the topographic humidity index in the study area, respectively, used for normalization. This represents the soil organic carbon content of the current pixel, in g / kg. The average soil organic carbon content in the study area is expressed in g / kg. The standard deviation of soil organic carbon content in the study area is expressed in g / kg. , and These are weighting coefficients, reflecting the contributions of slope, water accumulation, and carbon content to the source-sink pattern, respectively. Their values ​​are determined using the analytic hierarchy process (AHP) or expert scoring, and satisfy [the following criteria]. The normalization condition. In one embodiment of the present invention, the default value is taken as... , , The basis for this weighting is that slope is the primary factor controlling erosion, water accumulation determines the transport capacity of runoff, and spatial anomalies in carbon content reflect the cumulative effect of long-term erosion and deposition.

[0058] The physical meaning of this formula can be explained as follows: The first term reflects the influence of slope on the source-sink pattern. The greater the slope (above the critical slope), the stronger the erosion and the more likely it is to be a carbon export zone. The second term reflects the influence of water accumulation characteristics. The smaller the TWI (Temperature-Induced Width), the less likely water is to accumulate and the more dispersed the runoff, with erosion mainly consisting of dispersed surface erosion, and the more likely it is to be a carbon export zone. The larger the TWI, the more likely water is to accumulate and deposition may occur. The third term reflects the influence of spatial anomalies in carbon content. Areas with carbon content significantly higher than the mean may have organic carbon enrichment, while areas with carbon content significantly lower than the mean may have experienced carbon loss. The combined effect of these three terms comprehensively reflects the synergistic effect of topography, hydrology, and carbon distribution.

[0059] Then, the study area is divided into three functional zones based on the source-sink index. In one embodiment of the present invention, the division criteria are as follows: When SI is greater than 0.3, it is classified as an erosion carbon output zone. This area is located in the upper and middle part of the slope, with a large slope and dispersed surface runoff, and is the main source area for organic carbon to migrate downhill with eroded soil. When SI is between -0.3 and 0.3, it is classified as a transport transition zone. This area is located in the middle of the slope, where erosion and deposition processes occur alternately, and organic carbon is in a dynamic transport state. When SI is less than -0.3, it is classified as a deposition carbon enrichment zone. This area is located in the lower part of the slope or in the valley, with a large topographic humidity index, slow runoff velocity, and the organic carbon carried is enriched by sediment deposition. The threshold of 0.3 is selected based on: through field verification of multiple typical sample areas, when the absolute value of SI is greater than 0.3, its consistency rate with the actual erosion and deposition conditions exceeds 80%. This threshold can also be adjusted according to the actual situation of the study area, and is generally selected between 0.2 and 0.4.

[0060] Finally, the distribution map of carbon loss sources and sinks and the area statistics of each functional zone are output. Preferably, the differences in carbon content and carbon flux in each functional zone can be further calculated to analyze the spatial redistribution characteristics of organic carbon among the three functional zones. The carbon content in carbon-exporting zones is usually lower than the regional average, while the carbon content in carbon-enriched zones is usually higher than the regional average. This spatial heterogeneity reflects the lateral redistribution of organic carbon by hydraulic erosion.

[0061] Step S4, timing monitoring step.

[0062] The main purpose of this step is to differentiate the contribution of different seasonal erosion processes to carbon loss through comparative analysis of multi-temporal remote sensing images, thereby achieving spatiotemporal dynamic monitoring of carbon loss. In one embodiment of the present invention, this step specifically includes the following sub-steps.

[0063] First, multi-temporal hyperspectral remote sensing images are acquired before and after the snowmelt period and the rainy season. For high-altitude mountainous areas, the snowmelt period typically occurs from April to May, and the rainy season typically occurs from July to September. In one embodiment of this invention, images from the following four key temporal phases are selected: the pre-snowmelt period (late March to early April, when snow has not completely melted but bare soil is partially exposed), the post-snowmelt period (mid-to-late May, when snow has completely melted), the pre-rainy season period (late June to early July, when concentrated rainfall has not yet begun), and the post-rainy season period (mid-to-late September, when concentrated rainfall has basically ended). The multi-temporal images should ideally have similar solar altitude angles and vegetation phenological periods to minimize interference from non-erosive factors on carbon content changes.

[0064] Secondly, steps S1 to S3 are performed on the images of each time phase to obtain the spatial distribution map of soil carbon content, the coupling relationship between erosion and carbon loss, and the distribution of source and sink patterns for each time phase. In the time series analysis, the spectral inversion model can adopt a unified model of the same modeling period to ensure the comparability of inversion results at different times; alternatively, models can be established for each time phase, but systematic bias correction between models needs to be performed through common validation samples.

[0065] Then, the spatial distribution changes of soil carbon content in different time phases are compared. In one embodiment of the present invention, the calculation formulas for the changes in carbon content during the snowmelt period and the rainy season are as follows:

[0066] ,

[0067] ,

[0068] in: The value represents the change in carbon content during the snowmelt period, expressed in g / kg. A positive value indicates a decrease in carbon content during the snowmelt period, meaning that carbon loss has occurred. This represents the soil organic carbon content in the early stages of snowmelt, expressed in g / kg. This represents the soil organic carbon content during the later stages of snowmelt, expressed in g / kg. The value represents the change in carbon content during the rainy season, expressed in g / kg. A positive value indicates a decrease in carbon content during the rainy season, meaning that carbon loss has occurred. This represents the soil organic carbon content during the early rainy season, expressed in g / kg. This represents the soil organic carbon content during the later stages of the rainy season, expressed in g / kg.

[0069] Spatial distribution maps of carbon content changes can intuitively reflect the spatial patterns of carbon loss and carbon enrichment in the study area. Carbon loss during the snowmelt season is mainly caused by freeze-thaw cycles damaging soil structure and snowmelt runoff eroding topsoil, spatially manifested as significant carbon loss in the upper and middle parts of slopes and sunny slopes; carbon loss during the rainy season is mainly caused by rainfall splash and runoff erosion, spatially manifested as significant carbon loss on both sides of valleys and steep slopes.

[0070] Furthermore, the carbon loss due to snowmelt erosion and the carbon loss due to rainfall erosion are calculated separately. In one embodiment of the present invention, the formula for calculating the carbon loss per unit area is:

[0071] ,

[0072] ,

[0073] in: This represents carbon loss due to snowmelt erosion, expressed in g / m², indicating the net loss of organic carbon per unit area during the snowmelt period. This represents the amount of carbon lost due to rainfall erosion, expressed in g / m², indicating the net loss of organic carbon per unit area during the rainy season. This represents the change in carbon content during the snowmelt period, expressed in g / kg. This represents the change in carbon content during the rainy season, expressed in g / kg. This refers to the bulk density of the soil, expressed in g / cm³. Effective erosion depth, in cm; The carbon loss efficiency coefficient due to snowmelt erosion is dimensionless and ranges from 0.6 to 0.9. It reflects the efficiency of snowmelt runoff in carrying soil carbon. The basis for this coefficient is that snowmelt runoff has a relatively slow flow rate and a strong selective transport capacity for fine particulate organic carbon, but the total transport volume is limited by the amount of snowmelt water. The rainfall erosion carbon loss efficiency coefficient is a dimensionless coefficient ranging from 0.7 to 0.95. It reflects the efficiency of rainfall runoff in carrying soil carbon. This coefficient is usually slightly higher than that of snowmelt erosion because rainfall splash can more effectively disperse soil aggregates and release organic carbon.

[0074] Finally, the spatiotemporal dynamic monitoring results of carbon loss are output. The monitoring results include: spatial distribution maps of soil carbon content at different time phases, maps of carbon content changes during the snowmelt and rainy seasons, spatial distribution maps of carbon loss due to snowmelt erosion and rainfall erosion, and seasonal carbon flux statistics for each functional zone (carbon export zone, transition zone, and carbon enrichment zone). In one embodiment of the invention, the contribution rate of snowmelt erosion and rainfall erosion to the total annual carbon loss can also be calculated using the following formula:

[0075] ,

[0076] ,

[0077] in: The contribution rate of carbon loss due to snowmelt erosion represents the percentage of carbon loss during the snowmelt period relative to the total annual carbon loss. The contribution rate of precipitation erosion to carbon loss represents the percentage of carbon loss during the rainy season relative to the total annual carbon loss. For typical high-altitude mountainous areas, the contribution rate of snowmelt erosion is usually between 20% and 40%, while the contribution rate of precipitation erosion is usually between 60% and 80%, but the specific proportions vary depending on the regional climate characteristics and topographic conditions.

[0078] The above monitoring results can provide a scientific basis for carbon cycle research in high-altitude and cold mountainous areas and for assessing the carbon sink benefits of soil and water conservation.

[0079] See Figure 2 The present invention also provides a remote sensing image monitoring system for soil carbon loss in water erosion areas. This system corresponds to the above method embodiment and includes five core functional modules: a data acquisition module, a spectral inversion module 1, a coupling analysis module 2, a source-sink identification module 3, and a time series monitoring module 4.

[0080] The data acquisition module is used to acquire and manage various basic data required for monitoring, including hyperspectral remote sensing images of the study area, erosion intensity grading maps, digital elevation model data, and measured soil organic carbon content data from ground sampling points. In one embodiment of the present invention, the data acquisition module includes an image data interface, a vector data interface, and a field sampling data management submodule. The image data interface supports reading various hyperspectral satellite data formats, such as HDF5 format for GF-5 and TIFF format for HJ-1A / 1B; the vector data interface supports importing vector data such as erosion grading maps and land use maps; the field sampling data management submodule provides functions for inputting, storing, and querying sampling point coordinates, carbon content measurements, and auxiliary attribute information.

[0081] The spectral inversion module 1 is used to preprocess hyperspectral remote sensing images to obtain spectral reflectance data, extract spectral features in the visible-near-infrared bands, and construct a spectral inversion model using partial least squares regression and support vector machine algorithms to obtain soil carbon content inversion results. In one embodiment of the present invention, the spectral inversion module 1 includes a preprocessing submodule, a feature extraction submodule, and an inversion modeling submodule. The preprocessing submodule integrates radiometric correction, atmospheric correction, and geometric correction algorithms, and supports the selection of two atmospheric correction models, FLAASH and 6S. The feature extraction submodule implements spectral enhancement methods such as continuum removal, first-order differentiation, and discrete wavelet transform, and provides band sensitivity analysis and automatic feature band selection functions. The inversion modeling submodule integrates two modeling algorithms, PLSR and SVM, and provides functions such as parameter optimization, model training, cross-validation, and accuracy evaluation, supporting model saving, loading, and batch application. The spectral inversion module 1 executes the spectral inversion process described in step S1 of the method embodiment.

[0082] The coupling analysis module 2 is used to perform spatial overlay analysis between the soil carbon content inversion results and the erosion intensity classification map, calculate the statistical characteristic values ​​and carbon content difference values ​​of carbon content in areas with different erosion levels, and calculate the carbon loss rate based on the erosion intensity and carbon content difference values. In one embodiment of the present invention, the coupling analysis module 2 includes a spatial overlay submodule, a statistical analysis submodule, and a coupling relationship modeling submodule. The spatial overlay submodule realizes spatial registration of multiple layers and pixel-level overlay operations; the statistical analysis submodule calculates statistical indicators such as the mean, standard deviation, and coefficient of variation of carbon content according to erosion level partitions; the coupling relationship modeling submodule calculates the erosion-carbon loss coupling coefficient and carbon loss rate according to the formula in the method embodiment, and supports coupling relationship fitting of various function forms such as linear, exponential, and power functions. The coupling analysis module 2 executes the coupling analysis process described in step S2 of the method embodiment.

[0083] Source-sink identification module 3 is used to extract topographic parameters and construct a source-sink pattern discrimination index based on erosion and deposition characteristics, dividing the study area into erosion carbon output zone, migration transition zone, and deposition carbon enrichment zone. In one embodiment of the present invention, source-sink identification module 3 includes a topographic analysis submodule and a source-sink zoning submodule. The topographic analysis submodule extracts topographic parameters such as slope, aspect, runoff accumulation, topographic humidity index, and slope position index based on DEM data; the source-sink zoning submodule calculates the SI value of each pixel according to the source-sink index formula in the method embodiment, and divides functional areas according to a set threshold, outputting a source-sink pattern distribution map and a zoning statistical report. Source-sink identification module 3 executes the source-sink identification process described in step S3 of the method embodiment.

[0084] The time-series monitoring module 4 is used to acquire multi-temporal hyperspectral remote sensing images, compare the spatial distribution changes of soil carbon content in different time phases, and calculate the carbon loss due to snowmelt erosion and rainfall erosion, respectively. In one embodiment of the present invention, the time-series monitoring module 4 includes a multi-temporal image management submodule, a change detection submodule, and a carbon loss calculation submodule. The multi-temporal image management submodule supports the organization and management of hyperspectral images of different time phases according to time series, and provides image metadata browsing and time phase selection functions; the change detection submodule calls the spectral inversion module to perform carbon content inversion on each time phase image and calculates the carbon content change value between time phases; the carbon loss calculation submodule calculates the carbon loss during the snowmelt period and the rainy season according to the formula in the method embodiment, and generates a spatiotemporal dynamic monitoring report of carbon loss. The time-series monitoring module 4 executes the time-series monitoring process described in step S4 of the method embodiment.

[0085] Preferably, the system in this embodiment of the invention further includes a visualization module and a report generation module. The visualization module provides various visualization formats such as carbon content distribution maps, erosion-carbon loss relationship curves, source-sink pattern distribution maps, and time-series change animations, and supports layer overlay display and interactive querying; the report generation module automatically summarizes the analysis results of each module and generates a standardized monitoring report, which includes an overview of the study area, data sources, technical methods, monitoring results, and conclusions and recommendations.

[0086] The system provided in this invention can be deployed using a client-server architecture or a standalone architecture. In the client-server architecture, computationally intensive tasks such as spectral inversion and spatial analysis are executed on the server side, while the client provides user interaction and result display functions. In the standalone architecture, all functional modules are integrated into a desktop application, suitable for monitoring tasks in small to medium-sized study areas. System development can be based on Python combined with open-source libraries such as GDAL and Scikit-learn, or it can be based on IDL combined with the ENVI secondary development interface.

[0087] The embodiments of the present invention are not limited to the specific embodiments described above. Those skilled in the art can make various equivalent changes or substitutions based on the technical solutions of the present invention, and all such changes or substitutions should be included within the protection scope of the present invention.

Claims

1. A remote sensing image monitoring method for soil carbon loss in water erosion zones, characterized in that, The method includes the following steps: Step S1, Spectral Inversion Step: Acquire hyperspectral remote sensing images of the study area, preprocess the hyperspectral remote sensing images to obtain spectral reflectance data, extract spectral features in the visible-near-infrared bands, identify bare soil pixels, combine with measured soil organic carbon content data from ground sampling points, construct a spectral inversion model using partial least squares regression and support vector machine algorithms, perform areal estimation of soil organic carbon content in the surface layer of bare soil, and obtain soil carbon content inversion results; Step S2, Coupling Analysis Step: Obtain the erosion intensity classification map of the study area, perform spatial overlay analysis between the soil carbon content inversion results and the erosion intensity classification map, calculate the statistical characteristic value of carbon content and the carbon content difference value of different erosion level areas, calculate the carbon loss rate based on the erosion intensity and carbon content difference value, and establish the erosion-carbon loss coupling relationship. Step S3, source-sink identification step: Obtain digital elevation model data of the study area, extract topographic parameters such as topographic slope, topographic humidity index and runoff accumulation, construct source-sink pattern discrimination index in combination with erosion and deposition characteristics, divide the study area into erosion carbon output area, migration transition area and deposition carbon enrichment area, and identify the spatial redistribution pattern of organic carbon in the watershed. Step S4, temporal monitoring steps: acquire multi-temporal hyperspectral remote sensing images before and after the snowmelt period and before and after the rainy season, select images of four key time phases, execute steps S1 to S3 respectively, compare the spatial distribution changes of soil carbon content in different time phases, calculate the carbon loss due to snowmelt erosion and the carbon loss due to rainfall erosion respectively, and output the spatiotemporal dynamic monitoring results of carbon loss.

2. The remote sensing image monitoring method for soil carbon loss in water erosion zones according to claim 1, characterized in that, In step S1, the spectral features are extracted in the following bands: 480nm to 520nm, 580nm to 680nm, 1350nm to 1450nm, 1850nm to 1950nm, and 2150nm to 2250nm.

3. The remote sensing image monitoring method for soil carbon loss in water erosion zones according to claim 1, characterized in that, In step S1, bare soil pixels are identified using a threshold of NDVI less than 0.

2.

4. The remote sensing image monitoring method for soil carbon loss in water erosion areas according to claim 1, characterized in that, In step S2, the carbon content difference value is calculated based on the slightly eroded area, and the average carbon content difference between each level of erosion area and the slightly eroded area is calculated.

5. The remote sensing image monitoring method for soil carbon loss in water erosion areas according to claim 1, characterized in that, In step S2, the carbon loss rate is calculated by combining the erosion-carbon loss coupling coefficient, soil erosion modulus, soil bulk density, effective erosion depth, and soil organic carbon content.

6. The remote sensing image monitoring method for soil carbon loss in water erosion zones according to claim 1, characterized in that, In step S3, the source-sink pattern discrimination index comprehensively considers the difference between the terrain slope and the critical slope, the normalized value of the terrain humidity index, and the standardized deviation value of the soil organic carbon content.

7. The remote sensing image monitoring method for soil carbon loss in water erosion areas according to claim 6, characterized in that, In step S3, when the source-sink pattern discrimination index is greater than 0.3, it is classified as an erosion carbon output zone; when the source-sink pattern discrimination index is between -0.3 and 0.3, it is classified as a migration transition zone; and when the source-sink pattern discrimination index is less than -0.3, it is classified as a sedimentary carbon enrichment zone.

8. The remote sensing image monitoring method for soil carbon loss in water erosion zones according to claim 1, characterized in that, In step S4, images from four key time phases are selected: the pre-melting period, the post-melting period, the pre-rainy season period, and the post-rainy season period.

9. The remote sensing image monitoring method for soil carbon loss in water erosion zones according to claim 1, characterized in that, In step S4, the calculation of carbon loss due to snow melting erosion and carbon loss due to rainfall erosion is corrected by the carbon loss efficiency coefficient due to snow melting erosion and the carbon loss efficiency coefficient due to rainfall erosion, respectively. The carbon loss efficiency coefficient due to snow melting erosion ranges from 0.6 to 0.9, and the carbon loss efficiency coefficient due to rainfall erosion ranges from 0.7 to 0.

95.

10. A remote sensing image monitoring system for soil carbon loss in water erosion areas, used to implement the remote sensing image monitoring method for soil carbon loss in water erosion areas as described in any one of claims 1-9, characterized in that, include: The data acquisition module is used to acquire hyperspectral remote sensing images, erosion intensity classification maps, digital elevation model data, and measured data of soil organic carbon content at ground sampling points in the study area. The spectral inversion module is used to preprocess hyperspectral remote sensing images to obtain spectral reflectance data, extract spectral features in the visible-near-infrared bands, construct a spectral inversion model using partial least squares regression and support vector machine algorithms, and obtain soil carbon content inversion results. The coupling analysis module is used to perform spatial overlay analysis of soil carbon content inversion results and erosion intensity classification map, calculate the statistical characteristic value of carbon content and the carbon content difference value of different erosion level areas, and calculate the carbon loss rate based on erosion intensity and carbon content difference value. The source-sink identification module is used to extract topographic parameters and combine them with erosion and sedimentation characteristics to construct a source-sink pattern discrimination index, which divides the study area into an erosion carbon output zone, a migration transition zone, and a sedimentary carbon enrichment zone. The time-series monitoring module is used to acquire multi-temporal hyperspectral remote sensing images, compare the spatial distribution changes of soil carbon content in different time phases, and calculate the carbon loss due to snowmelt erosion and rainfall erosion, respectively.

Images