GRACE satellite-based drought region identification method

Abstract

The invention discloses an arid region identification method based on a GRACE satellite, and the method comprises the steps: obtaining satellite remote sensing original data, carrying out the median filtering of the satellite remote sensing original data, and obtaining noise reduction data; complementing missing data for the noise reduction data to obtain complemented data; calculating a normalized vegetation index and a temperature vegetation drought index according to the complementation data to obtain a drought index; time sequence features are extracted according to the drought indexes, principal component analysis is carried out on the time sequence features, and comprehensive indexes are obtained; spatial clustering analysis is carried out on the comprehensive indexes, and target areas with similar drought characteristics are divided; calculating a topographic relief degree and a land utilization type proportion based on the target area in combination with the land utilization data and the topographic data; selecting a target area of which the topographic relief degree is greater than a preset threshold value, and performing topographic correction on the drought index of the selected target area to obtain a corrected drought index; and judging the drought grade of the target area according to the corrected drought index and the land utilization type ratio.

Description

Technical Field

[0001] The present invention belongs to the technical field of arid region identification, and particularly relates to a method for identifying arid regions based on GRACE satellites. Background Art

[0002] Satellite remote sensing technology plays an important role in arid region identification, but faces challenges in data quality and processing methods. First, the original data obtained by satellites often have problems such as noise, outliers, and data missing, which need to be preprocessed to improve the accuracy and reliability of the data. However, different preprocessing methods may lead to data changes and distortions, affecting the results of subsequent analysis. Second, the selection and calculation method of drought indices are also a key issue. Different drought indices reflect different characteristics of drought, such as soil moisture, vegetation conditions, etc., and appropriate indices need to be selected according to the characteristics of the study area. However, the calculation of drought indices often involves multiple parameters and assumptions, and different parameter combinations and threshold settings may lead to differences in identification results. In addition, the identification of arid regions also needs to consider the issues of time scale and space scale. Drought is a continuous process, and it is necessary to analyze the drought change trend within a certain time range; at the same time, drought is regional and uneven in space, and it is necessary to comprehensively consider the environmental conditions and land use status of different regions. Therefore, how to select appropriate spatio-temporal scales for arid region identification is a challenging problem, which requires balancing the contradiction between the spatio-temporal resolution of data and the computational efficiency. Summary of the Invention

[0003] To solve the above technical problems, the present invention proposes a method for identifying arid regions based on GRACE satellites, which improves the accuracy and reliability of drought monitoring.

[0004] To achieve the above object, the present invention provides a method for identifying arid regions based on GRACE satellites, including: acquiring satellite remote sensing original data, performing median filtering processing on the satellite remote sensing original data to obtain noise-reduced data;

[0005] Completing the missing data of the noise-reduced data to obtain completed data;

[0006] Calculating the normalized difference vegetation index and temperature vegetation drought index according to the completed data to obtain drought indicators;

[0007] Extracting time series features according to the drought indicators, and performing principal component analysis on the time series features to obtain comprehensive indicators;

[0008] Performing spatial clustering analysis on the comprehensive indicators to divide target regions with similar drought characteristics;

[0009] Based on the target area, combine land use data and terrain data to calculate the terrain undulation degree and the proportion of land use types;

[0010] Select the target areas where the terrain undulation degree is greater than the preset threshold, and perform terrain correction on the drought index of the selected target areas to obtain the corrected drought index;

[0011] Determine the drought level of the target area according to the corrected drought index and the proportion of land use types.

[0012] Technical effects of the present invention: The present invention discloses a method for identifying drought areas based on GRACE satellites. First, noise reduction and interpolation processing are performed on the original data, and then the normalized difference vegetation index and the temperature vegetation drought index are calculated. Time series features are extracted and principal component analysis is performed. Target areas with similar drought characteristics are divided through spatial clustering analysis, and correction is performed in combination with land use and terrain data. Finally, a decision tree algorithm is used to determine the drought level. The present invention solves the deficiencies of traditional drought monitoring methods in data processing, feature extraction, and regional division, and improves the accuracy and reliability of drought monitoring. The present invention can effectively identify and evaluate drought conditions under different terrain conditions, providing an important basis for water resource management and agricultural production decision-making. Description of the Drawings

[0013] The drawings forming a part of this application are used to provide a further understanding of this application. The schematic embodiments of this application and their descriptions are used to explain this application and do not constitute an improper limitation to this application. In the drawings:

[0014] Figure 1 It is a flow schematic diagram of the method for identifying drought areas based on GRACE satellites in the embodiment of the present invention. Detailed Embodiments

[0015] It should be noted that, without conflict, the embodiments in this application and the features in the embodiments can be combined with each other. The following will refer to the drawings and combine the embodiments to detail this application.

[0016] It should be noted that the steps shown in the flowchart of the drawings can be executed in a computer system such as a set of computer executable instructions, and although the logical order is shown in the flowchart, in some cases, the steps shown or described can be executed in a different order than here.

[0017] As Figure 1 shown, the method for identifying drought areas based on GRACE satellites provided in this embodiment includes:

[0018] Obtain the original satellite remote sensing data, perform median filtering on the original satellite remote sensing data to obtain noise-reduced data;

[0019] Complement the missing data in the noise-reduced data to obtain the completed data;

[0020] Calculate the normalized difference vegetation index and temperature vegetation drought index based on the completed data to obtain drought indicators;

[0021] Extract the time series features according to the drought indicators, and perform principal component analysis on the time series features to obtain comprehensive indicators;

[0022] Perform spatial clustering analysis on the comprehensive indicators to divide the target areas with similar drought characteristics;

[0023] Based on the target areas, combine the land use data and topographic data to calculate the topographic undulation degree and the proportion of land use types;

[0024] Select the target areas where the topographic undulation degree is greater than the preset threshold, and perform topographic correction on the drought indicators of the selected target areas to obtain the corrected drought indicators;

[0025] Determine the drought level of the target areas according to the corrected drought indicators and the proportion of land use types.

[0026] Further, obtaining the noise-reduced data includes:

[0027] Determine the window size parameter of the median filtering algorithm according to the remote sensing data type;

[0028] Process the original satellite remote sensing data with the median filtering algorithm with the determined parameters to obtain the filtered data;

[0029] Judge whether the filtered data meets the preset noise reduction standard. If it does not meet the noise reduction standard, adjust the filtering window size parameter and reprocess to obtain the noise-reduced data; if it meets the noise reduction standard, directly obtain the noise-reduced data.

[0030] Further, obtaining the completed data includes:

[0031] Analyze the spatial distribution characteristics of the data missing areas according to the noise-reduced data;

[0032] Determine the parameter range of spatio-temporal interpolation according to the spatial distribution characteristics of the missing areas;

[0033] Perform spatial interpolation processing on the missing areas. If there are outliers in the interpolation results, correct the outliers;

[0034] Judge whether the completed data meets the preset integrity standard. If it does not meet the integrity standard, adjust the interpolation parameters and reprocess the missing areas to obtain the completed data.

[0035] Specifically, determining whether the complemented data meets the preset integrity criteria is a crucial step in ensuring data quality. Integrity criteria may include indicators such as data coverage rate and spatial continuity. For example, it may be required that in any 50×50 pixel area of the complemented coastline image, the proportion of valid data is not less than 95%. If it is found that some areas still do not meet this standard, it is necessary to return and adjust the interpolation parameters and reprocess these areas. Finally, the complemented remote sensing data is stored as a new complete dataset. This dataset not only fills in the missing parts of the original data but also undergoes outlier correction and integrity verification, providing high-quality basic data for subsequent analysis and applications. For example, this complete remote sensing image of the coastline can be used in multiple fields such as coastline change monitoring, marine resource assessment, or climate change research, providing important basis for relevant decision-making.

[0036] Furthermore, obtaining the drought index includes:

[0037] Based on the complemented data, extract the reflectance values of the near-infrared band and the red band, and calculate the vegetation index;

[0038] Extract the radiance value of the thermal infrared band and retrieve the land surface temperature;

[0039] Based on the vegetation index and the land surface temperature, calculate the drought index;

[0040] Determine whether the drought index is within the preset range. If it exceeds the range, adjust the formula parameters and recalculate; if the drought index meets the preset range, directly obtain the drought index.

[0041] Specifically, after obtaining the complementary data, first extract the reflectance values in the near-infrared and red light bands. These two bands are particularly important for vegetation monitoring because vegetation exhibits unique reflection characteristics in these bands. For example, healthy vegetation has a high reflectance in the near-infrared band and a low reflectance in the red light band. Using this property, the Normalized Difference Vegetation Index (NDVI) can be calculated. The NDVI value usually ranges from -1 to 1, where a higher value indicates a higher vegetation coverage, and a lower value may indicate bare soil or water bodies. Extract the radiance values from the thermal infrared band. The thermal infrared band can reflect the thermal radiation characteristics of the earth's surface. Through the Planck formula, the radiance can be converted into the surface temperature. This step is very important for studying the thermal environment of the research area. For example, in the study of the urban heat island effect, it may be found that the surface temperature in the urban central area is significantly higher than that in the surrounding suburbs. Combine the vegetation index and the surface temperature to calculate the Temperature-Vegetation Dryness Index (TVDI). TVDI is an important indicator for evaluating the drought degree of a region, comprehensively considering the vegetation coverage and temperature factors. For example, under the same temperature conditions, areas with lower vegetation coverage may be more vulnerable to drought. Carefully select the parameters of the TVDI formula to ensure that the calculation results can accurately reflect the actual situation of the research area. The calculated drought index needs to be verified to ensure that it is within the preset range.

[0042] Furthermore, dividing the target areas with similar drought characteristics includes:

[0043] Standardize the comprehensive indicators to eliminate the dimensional differences between different indicators;

[0044] Conduct a spatial autocorrelation analysis on the standardized comprehensive indicators, calculate the global Moran's I index, and judge the spatial distribution pattern of the comprehensive indicators;

[0045] According to the results of the spatial autocorrelation analysis, select the spatial weight matrix and construct a spatial relationship model;

[0046] Conduct a spatially constrained clustering analysis on the comprehensive indicators, construct a spatial adjacency relationship, and divide the clusters based on the dendrogram cutting principle to obtain the clustering results;

[0047] Check the spatial continuity of the clustering results. If there are spatial fragments, merge the fragments to ensure the spatial continuity of the clustering areas;

[0048] Calculate the mean and standard deviation of the comprehensive indicators within each clustering area, and rank the clustering areas according to the mean value to obtain the division results of the target areas with similar drought characteristics.

[0049] Specifically, after obtaining the comprehensive index data from the database, the first step is to perform standardization, which is crucial. Suppose there are two indicators: precipitation and temperature. The unit of precipitation is millimeters, and the unit of temperature is degrees Celsius. Directly comparing these two indicators will result in deviations due to different units. Through standardization, the two indicators can be converted to the same scale, such as converting to standard scores (Z-scores). After such processing, the values between different indicators become comparable. Next, spatial autocorrelation analysis is carried out to calculate the global Moran's I. This step helps to understand the spatial distribution pattern of drought characteristics. Suppose the global Moran's I of a certain region is 0.85, which indicates that the drought characteristics in this region show strong spatial clustering. Areas with similar drought degrees tend to cluster together rather than being randomly distributed. Based on the results of spatial autocorrelation analysis, an appropriate spatial weight matrix is selected. If the analysis shows that the drought characteristics are similar in a large area, an inverse distance weight matrix based on distance may be chosen. This matrix assigns higher weights to regions closer in distance, thus considering spatial proximity relationships in subsequent clustering analysis. The application of the SKATER algorithm is the core of this process. The algorithm represents spatial adjacency relationships by constructing a minimum spanning tree and then performs clustering based on the principle of branch cutting. Suppose there is a province composed of 100 county-level administrative regions. The SKATER algorithm will first construct a minimum spanning tree connecting all county-level units, and then divide these 100 counties into several regions with similar drought characteristics according to the similarity of comprehensive indicators and spatial adjacency relationships.

[0050] Furthermore, calculating the terrain undulation degree and the proportion of land use types includes:

[0051] Obtain the land use data and terrain data of the target area, preprocess the land use data and terrain data to obtain the preprocessed data;

[0052] Perform spatial interpolation on the preprocessed data to generate continuous land use type distribution maps and terrain undulation degree maps;

[0053] According to the terrain undulation degree map, divide the terrain grades and calculate the area proportion of each grade;

[0054] Obtain the land use type distribution map and count the area proportion of each land use type.

[0055] Furthermore, obtaining the corrected drought index includes:

[0056] According to the terrain undulation degree map, extract the terrain undulation degree value of each grid cell, compare it with a pre-set terrain undulation degree threshold, and determine whether the grid cell needs to perform terrain correction on the drought index;

[0057] According to the grid cells that need to be corrected as judged, obtain the corresponding original drought index values, and use a terrain correction model to correct the drought index to obtain the corrected drought index.

[0058] Specifically, use geographic information system software such as ArcGIS to calculate the slope and aspect based on elevation data. The slope reflects the steepness of the terrain, and the aspect indicates the orientation of the slope surface. These two parameters have an important impact on the degree of drought. When constructing a terrain undulation model, a neighborhood analysis method can be used. Taking a moving window of 500 meters × 500 meters as an example, calculate the difference between the maximum elevation and the minimum elevation within the window as the terrain undulation. This method can better reflect the undulation changes of the local terrain. The obtained terrain undulation distribution map can visually display the complexity of the terrain within the region. Terrain correction is a key step in improving the accuracy of drought indices. Taking the normalized difference vegetation index (NDVI) as an example, when the terrain undulation exceeds a preset threshold (such as 100 meters), correction is required. The correction model can consider the impact of slope on vegetation growth. For example, in steep slope areas, vegetation growth is poor and the NDVI value is low, so it needs to be appropriately increased. At the same time, the aspect also affects vegetation growth. For example, in the Northern Hemisphere, the south slope has sufficient sunlight and better vegetation growth, and the NDVI value may be overestimated and needs to be appropriately reduced. Spatial interpolation is an important method for generating continuous distribution maps. Kriging interpolation takes into account spatial autocorrelation and is suitable for regions with large terrain changes. For example, in the mountain-plain transition zone, ordinary Kriging is used for interpolation at a grid spacing of 500 meters to obtain a continuous corrected drought index distribution map.

[0059] Furthermore, determining the drought level of the target area includes:

[0060] Based on the corrected drought index and the proportion of land use types, perform spatial registration and resampling processing to obtain the same spatial resolution and geographic coordinates;

[0061] Using the same spatial resolution and geographic coordinates, calculate the information gain ratio of each feature, and select the feature with the largest information gain ratio as the root node of the decision tree to determine the optimal splitting feature;

[0062] According to the determined optimal splitting feature, recursively construct a decision tree model. Calculate the information gain ratio of the child nodes at each non-leaf node, and select the optimal splitting feature and splitting point until the stopping condition is reached;

[0063] At the leaf nodes of the decision tree, set the drought level determination threshold according to the class distribution of the samples, and divide the leaf nodes into different levels to form drought level determination rules;

[0064] Perform pruning on the decision tree model. By comparing the error rates before and after pruning, remove the overfitting branches to obtain the pruned decision tree model;

[0065] Input the drought index and land use type of the target area into the pruned decision tree model, determine the drought level for each grid cell, obtain the grid-based drought level distribution result, and generate the drought level distribution map of the target area.

[0066] Specifically, taking a plain rice-growing area as an example, MODIS Normalized Difference Vegetation Index (NDVI) data with a 1-kilometer resolution and land use type data with a 30-meter resolution were obtained. To achieve the matching of the two types of data, the nearest neighbor resampling method was used to resample the land use data to a 1-kilometer resolution and unified to the WGS84 coordinate system. This process ensured the spatial consistency of subsequent analyses. Feature extraction is the basis for constructing a decision tree model. Vegetation health status indicators such as the annual change rate of NDVI and the mean NDVI during the growing season can be extracted from the NDVI data. Features such as the proportion of farmland and the proportion of forest land can be extracted from the land use data. For example, calculate the proportion of paddy field area in each 1-square-kilometer grid to reflect the regional water resource utilization situation. These features together constitute a multi-dimensional vector describing the regional drought vulnerability. Feature selection algorithms aim to improve the efficiency and accuracy of the model. The calculation of the information gain ratio can reveal the contribution degree of each feature to the drought level determination. Suppose the calculation result shows that the information gain ratio of the annual change rate of NDVI is the highest, reaching 0.8, then it is selected as the root node of the decision tree. This selection reflects the importance of vegetation dynamic changes in drought monitoring. The decision tree construction process reflects the hierarchy of drought assessment. Taking the C5 algorithm as an example, the possible splitting paths are as follows: First, divide according to the annual change rate of NDVI into normal years (>-5%) and abnormal years (≤-5%); for abnormal years, further divide according to the proportion of paddy fields into high-risk (>60%) and medium-low-risk (≤60%) areas. This hierarchical division reflects the comprehensive impact of vegetation conditions and land use types on drought sensitivity. The drought level determination of leaf nodes reflects the flexibility of classification. For example, a leaf node with an annual change rate of NDVI < -10% and a paddy field proportion > 80% may be determined as a severe drought risk area. This determination takes into account the dual factors of severe vegetation degradation and high dependence on irrigation. Decision tree pruning is an important means to improve the generalization ability of the model. Suppose the determination accuracy of a certain branch is only 55%, lower than the overall accuracy of 70%, then it may be pruned. This process helps to reduce the overfitting of the model to the training data and improve its performance on new data. The finally generated drought level distribution map visually shows the spatial pattern of regional drought risk. For example, the map may show a low-risk zone along the riverbank and a high-risk area far from the water source. This visualization result provides an important basis for water resource management and agricultural production planning, and helps to formulate targeted drought resistance measures.

[0067] The present invention discloses a method for identifying drought regions based on GRACE satellites. First, noise reduction and interpolation processing are performed on the original data, then the normalized difference vegetation index and the temperature vegetation drought index are calculated, time series features are extracted and principal component analysis is carried out. Target regions with similar drought characteristics are divided through spatial clustering analysis, and correction is performed in combination with land use and terrain data. Finally, a decision tree algorithm is used to determine the drought level. The present invention solves the deficiencies of traditional drought monitoring methods in data processing, feature extraction, and regional division, and improves the accuracy and reliability of drought monitoring. The present invention can effectively identify and evaluate drought conditions under different terrain conditions, providing an important basis for water resource management and agricultural production decision-making.

[0068] The above is only a preferred specific embodiment of the present application, but the protection scope of the present application is not limited thereto. Any changes or substitutions that can be easily thought of by those skilled in the art within the technical scope disclosed in the present application should be covered by the protection scope of the present application. Therefore, the protection scope of the present application should be subject to the protection scope of the claims.

Claims

1. A drought area identification method based on GRACE satellite, characterized in that: include: Acquire raw satellite remote sensing data, perform median filtering on the raw satellite remote sensing data, and obtain noise-reduced data; Completing missing data of the noise reduction data to obtain completed data; Calculating the normalized vegetation index and the temperature vegetation drought index according to the supplemented data to obtain a drought index; Extracting time series characteristics according to the drought index, performing principal component analysis on the time series characteristics, and obtaining a comprehensive index; Performing spatial cluster analysis on the comprehensive index to divide target areas with similar drought characteristics; Calculate the terrain relief and the proportion of land use types based on the target area in combination with the land use data and the terrain data; Selecting a target area where the terrain relief is greater than a preset threshold, performing terrain correction on the drought index of the selected target area to obtain a corrected drought index; The drought level of the target area is determined according to the corrected drought index and the land use type ratio.

2. The drought area identification method based on GRACE satellite according to claim 1, characterized in that: Obtaining noise reduction data includes: Determine the window size parameter of the median filter algorithm according to the remote sensing data type; Processing the satellite remote sensing raw data using a median filtering algorithm with determined parameters to obtain filtered data; It is determined whether the filtered data meets the preset noise reduction standard. If it does not meet the noise reduction standard, the filter window size parameter is adjusted and reprocessed to obtain the noise reduction data; if it meets the noise reduction mark, the noise reduction data is directly obtained.

3. The drought area identification method based on GRACE satellite according to claim 1, characterized in that: The completed data includes: Analyzing spatial distribution characteristics of data missing areas based on the noise reduction data; Determine the parameter range of spatiotemporal interpolation based on the spatial distribution characteristics of the missing area; Perform spatial interpolation on the missing area, and correct the outliers if there are any in the interpolation results; It is determined whether the completed data meets the preset integrity standard. If it does not meet the integrity standard, the interpolation parameters are adjusted to reprocess the missing area to obtain the completed data.

4. The drought area identification method based on GRACE satellite according to claim 1, characterized in that: Obtaining drought indicators includes: Based on the supplemented data, the reflectance values ​​of the near infrared band and the red light band are extracted to calculate the vegetation index; Extract the radiation brightness value of the thermal infrared band and invert the surface temperature; Calculating a drought index based on the vegetation index and the land surface temperature; Determine whether the drought index is within the preset range. If it is beyond the range, adjust the formula parameters and recalculate; if the drought index meets the preset range, directly obtain the drought index.

5. The drought area identification method based on GRACE satellite according to claim 1, characterized in that: Target areas with similar drought characteristics include: Standardizing the comprehensive indicators to eliminate the dimensional differences between different indicators; The spatial autocorrelation analysis of the standardized comprehensive index was carried out, and the global Moran index was calculated to determine the spatial distribution pattern of the comprehensive index; According to the results of spatial autocorrelation analysis, the spatial weight matrix is ​​selected to construct the spatial relationship model; Perform spatial restricted cluster analysis on the comprehensive indicators, construct spatial adjacency relationships, and divide clusters based on the tree branch cutting principle to obtain clustering results; Performing a spatial continuity check on the clustering results, and if there are spatial fragments, merging the fragments to ensure the spatial continuity of the clustering area; The mean and standard deviation of the comprehensive indicators in each cluster area were calculated, and the cluster areas were ranked according to the mean value to obtain the division results of target areas with similar drought characteristics.

6. The drought area identification method based on GRACE satellite according to claim 1, characterized in that: Calculation of terrain relief and land use type ratio includes: Acquire land use data and terrain data of the target area, preprocess the land use data and terrain data, and obtain preprocessed data; Performing spatial difference processing on the preprocessed data to generate a continuous land use type distribution map and a terrain relief map; According to the terrain relief map, the terrain is divided into grades and the area proportion of each grade is calculated; Obtain a land use type distribution map and count the area proportion of each land use type.

7. The drought area identification method based on GRACE satellite according to claim 1, characterized in that: Corrected drought indicators include: According to the terrain relief map, the terrain relief value of each grid cell is extracted and compared with the pre-set terrain relief threshold to determine whether the grid cell needs terrain correction of drought indicators; According to the grid cells that need to be corrected, the corresponding original drought index value is obtained, and the drought index is corrected using the terrain correction model to obtain the corrected drought index.

8. The drought area identification method based on GRACE satellite according to claim 1, characterized in that: Determining the drought level of the target area includes: According to the calibrated drought index and land use type ratio, spatial registration and resampling are performed to obtain the same spatial resolution and geographic coordinates; Using the same spatial resolution and geographic coordinates, calculate the information gain ratio of each feature, select the feature with the largest information gain ratio as the root node of the decision tree, and determine the optimal split feature; According to the determined optimal splitting features, a decision tree model is recursively constructed, the information gain ratio of the child nodes is calculated at each non-leaf node, and the optimal splitting features and splitting points are selected until the stopping condition is reached; At the leaf nodes of the decision tree, the drought level determination threshold is set according to the category distribution of the samples, and the leaf nodes are divided into different levels to form drought level determination rules; Pruning the decision tree model, removing overfitted branches by comparing the error rates before and after pruning, and obtaining a pruned decision tree model; The drought index and land use type of the target area are input into the pruned decision tree model, the drought level is determined for each grid unit, a gridded drought level distribution result is obtained, and a drought level distribution map of the target area is generated.

Images