A method for constructing agricultural drought index considering spatial heterogeneity

Abstract

The present invention provides a method for constructing an agricultural drought index considering spatial heterogeneity, including constructing a generalized agricultural drought index GSMSDI, and also including the following steps: determining the agricultural drought index SMSDI corresponding to different thresholds; obtaining MODIS13A3 products, extracting monthly scale normalized vegetation index data; obtaining MOD12Q1 products, extracting land cover data sets of each year; using an optimal threshold extraction method to extract the optimal threshold for defining the drought index pixel by pixel; obtaining the optimal threshold for defining the agricultural drought index under different underlying surfaces. The present invention proposes a method for constructing an agricultural drought index considering spatial heterogeneity, based on time series soil moisture, constructing a drought index with the help of vegetation index data, so that the index has a significant agricultural drought indication effect.

Description

Technical Field

[0001] The invention relates to the technical field of drought monitoring, and in particular to a method for constructing an agricultural drought index taking spatial heterogeneity into consideration. Background Art

[0002] Drought index is an index indicating drought status. There are many drought indices. In the monitoring application of agricultural drought, the index can usually be divided into precipitation-based index, soil moisture or soil water balance-based index and vegetation growth status-based index. These indices focus on a certain stage or part of the agricultural drought process (precipitation-soil moisture-vegetation growth status-crop yield).

[0003] Among the above three types of agricultural drought indices, the precipitation-based index usually defines the degree of drought by precipitation. The indicated drought conditions cannot be effectively transmitted and ultimately affect the final crop yield. This is because, on the one hand, precipitation will be lost through runoff and other means; on the other hand, affected by soil texture, vegetation root distribution, etc., different types of vegetation have different sensitivities to precipitation responses. The drought index defined based on vegetation growth conditions effectively reflects the impact of drought, but from the perspective of water management, when vegetation has been affected by drought, it has limited significance for guiding water management and alleviating agricultural drought. Therefore, the index based on soil moisture is usually more applicable. First, soil moisture can effectively affect the growth of vegetation and can be used to predict crop yields; second, by monitoring soil moisture, it can effectively warn of the occurrence of drought, which is crucial for water management. Given that long-term soil moisture data can be widely obtained, the development of drought indices based on soil moisture time series statistical data can be more widely used.

[0004] At present, there is no unified definition of drought occurrence in the index based on soil moisture development, that is, the threshold is uncertain when defining the index. The commonly used index determination method is usually based on historical statistical data, and the index is constructed with a certain quantile as the threshold (such as the 10th quantile of historical time series data). The index constructed by this method is comparable in space, but considering that the drought indicated by the drought index is ultimately manifested as crop yield, therefore, it is only a good agricultural drought index when the index can well predict crop yield or vegetation productivity. This shows that the definition of drought index is still heterogeneous in space due to differences in vegetation response to drought. Therefore, taking these differences into full consideration and extracting a more reasonable drought index is of great significance for drought early warning, crop yield prediction, etc.

[0005] The invention patent application with publication number CN113919146A discloses a method for constructing an agricultural drought index based on soil available water, which uses the WDt index to construct the agricultural drought index. Based on the root layer soil water content simulated by the crop growth model (WOFOST), a modified soil water deficit index (MSWDI) suitable for agricultural drought monitoring was constructed. Although this method also takes into account vegetation type (such as root depth, etc.), it requires simulation parameters and soil texture parameters, which are difficult to obtain using remote sensing technology and have large uncertainties on a macro scale. .

[0006] The invention patent application with publication number CN113095621A discloses an agricultural drought detection method based on the time lag of soil moisture to meteorological conditions, and develops the CADIi drought index, which simultaneously considers three parameters: precipitation, evapotranspiration and soil moisture. The final index has a better comparison effect with the SPEI, but does not consider the differences in vegetation response to drought, that is, it does not consider spatial heterogeneity like the present invention. Summary of the invention

[0007] In order to solve the above technical problems, the present invention proposes a method for constructing an agricultural drought index taking into account spatial heterogeneity. Based on time series soil moisture and with the help of vegetation index data, a drought index is constructed so that the index has a significant agricultural drought indication effect.

[0008] The present invention provides a method for constructing an agricultural drought index considering spatial heterogeneity, including constructing a generalized agricultural drought index GSMSDI, and further comprising the following steps:

[0009] Step 1: Determine the agricultural drought index SMSDI corresponding to different thresholds;

[0010] Step 2: Obtain MODIS13A3 products and extract monthly normalized difference vegetation index data;

[0011] Step 3: Obtain MOD12Q1 products and extract land cover datasets for each year;

[0012] Step 4: Use the optimal threshold extraction method to extract the optimal threshold for defining the drought index pixel by pixel;

[0013] Step 5: Obtain the optimal threshold for defining the agricultural drought index under different underlying surfaces.

[0014] Preferably, the generalized agricultural drought index GSMSDI is constructed based on the time series soil moisture, and the generalized agricultural drought index GSMSDI is defined as

[0015]

[0016] Among them, P smis the percentile value of soil moisture in the time series data, T0 is the percentile threshold, and the generalized meaning is that the percentile T0 has not been determined at this time.

[0017] In any of the above schemes, it is preferred that when the soil moisture percentile value is lower than T0, GSMSDI<0, when the soil moisture percentile value is higher than T0, GSMSDI>0, otherwise, GSMSDI=0.

[0018] In any of the above schemes, it is preferred that the GSMSDI index ranges from -1 to 1, the closer to 1, the wetter it is, and the closer to -1, the drier it is.

[0019] In any of the above schemes, preferably, step 1 includes selecting T0=5, 6, 7,…, 50th for the generalized agricultural drought index GSMSDI, taking values ​​every step of 1, and determining the agricultural drought index SMSDI corresponding to different thresholds, thereby generating a total of 46 time series.

[0020] In any of the above schemes, preferably, step 2 includes calculating the NDVI anomaly index NDVIA using NDVI data, and the formula of the NDVI anomaly index NDVIA is:

[0021]

[0022] Among them, NDVIA i,j is the NDVI anomaly index of the jth month of the i-th year, i is the year sequence, j is the jth period of the growing season, the growing season is from April to September, and the NDVI ave,j is the multi-year NDVI mean of the jth period.

[0023] In any of the above solutions, preferably, step 2 further includes obtaining the annual NDVIA by calculating the monthly average value after calculating the monthly NDVIA.

[0024] In any of the above schemes, it is preferred that the optimal threshold extraction method is to calculate the correlation between multi-year time series SMSDI data and NDVIA data under different threshold conditions, and the threshold corresponding to when the correlation between the two reaches the maximum value is used to define the optimal threshold of the drought index.

[0025] In any of the above solutions, preferably, the calculation formula of the correlation r is:

[0026]

[0027] Among them, SMSDI i SMSDI, NDVIA for the i-th year i is the NDVIA for the i-th year; is the mean value of SMSDI, are the mean values ​​of NDVIA; r max is the maximum correlation coefficient, and n is the total number of years.

[0028] In any of the above schemes, preferably, the maximum correlation coefficient r max The calculation formula is

[0029] r max =max{r T0＝5th ,r T0＝6th ,r T0＝7th ,…,r T0＝50th}

[0030] Among them, r T0＝5th ,r T0＝6th ,r T0＝7th ,…,r T0＝50th is the correlation coefficient between SMSDI and NDVIA calculated when T0 takes the 5th, 6th, 7th, …, 50th percentile values ​​respectively, and max is the function of taking the maximum value in the set.

[0031] In any of the above schemes, preferably, step 5 includes the following sub-steps:

[0032] Step 51: Determine the optimal threshold for each pixel;

[0033] Step 52: combining land cover data, calculating the mean value of the optimal threshold of all pixels under the specific underlying surface;

[0034] Step 53: Determine the optimal threshold value defined for the drought index of a specific underlying surface using the mean value.

[0035] In any of the above solutions, preferably, the step 51 is to determine the optimal threshold pixel by pixel using the optimal threshold extraction method.

[0036] MODIS13A3 products refer to vegetation index products provided by the Terra / Aqua medium-resolution imaging spectroradiometer, including the Normalized Difference Vegetation Index (NDVI), with a spatial resolution of 1 km and a temporal resolution of months. Other satellite or reanalysis vegetation index products can also be used, as long as the corresponding vegetation index can reflect the vegetation growth status.

[0037] MOD12Q1 products refer to land cover products provided by the Terra / Aqua medium-resolution imaging spectroradiometer, which include land cover datasets with a spatial resolution of 1 km and one land cover classification image per year. Other satellite or reanalysis land cover datasets can also be used to obtain land cover types.

[0038] NDVI, Normalized Difference Vegetation Index, refers to the normalized vegetation index on a monthly scale. BRIEF DESCRIPTION OF THE DRAWINGS

[0039] Figure 1 The present invention is a flowchart of a preferred embodiment of a method for constructing an agricultural drought index taking spatial heterogeneity into consideration.

[0040] Figure 2 This is a schematic diagram of SMSDI changes corresponding to different thresholds of soil moisture extracted from a certain time series according to a free embodiment of the method for constructing an agricultural drought index taking into account spatial heterogeneity of the present invention.

[0041] Figure 3 This is a schematic diagram of the variation of the correlation coefficient between NDVIA and SMSDI at a certain pixel with different SMSDI thresholds according to a preferred embodiment of the method for constructing an agricultural drought index taking into account spatial heterogeneity of the present invention.

[0042] Figure 4 A statistical schematic diagram of the optimal threshold pixel number at a regional scale according to a preferred embodiment of a method for constructing an agricultural drought index taking spatial heterogeneity into consideration according to the present invention.

[0043] Figure 5 This is a schematic diagram comparing the SPEI and SMSDI change curves of three stations according to a preferred embodiment of the method for constructing an agricultural drought index taking into account spatial heterogeneity of the present invention. DETAILED DESCRIPTION

[0044] The present invention is further described below in conjunction with the accompanying drawings and specific embodiments.

[0045] Embodiment 1

[0046] like Figure 1 As shown, step 100 is executed to construct a generalized agricultural drought index GSMSDI based on the time series soil moisture. The generalized agricultural drought index GSMSDI is defined as

[0047]

[0048] Among them, P sm is the percentile value of soil moisture in the time series data, T0 is the percentile threshold, and the generalized meaning is that the percentile T0 has not yet been determined. When the soil moisture percentile value is lower than T0, GSMSDI < 0, when the soil moisture percentile value is higher than T0, GSMSDI > 0, otherwise, GSMSDI = 0. The GSMSDI index ranges from -1 to 1. The closer it is to 1, the wetter it is, and the closer it is to -1, the drier it is.

[0049] Execute step 110 to determine the agricultural drought index SMSDI corresponding to different thresholds. For the generalized agricultural drought index GSMSDI, select T0=5, 6, 7, ..., 50th, take values ​​every step 1, and determine the agricultural drought index SMSDI corresponding to different thresholds, thereby generating a total of 46 time series.

[0050] Execute step 120 to obtain MODIS13A3 products and extract monthly normalized difference vegetation index data. NDVI anomaly index NDVIA is calculated using the following formula using NDVI data. The formula for NDVI anomaly index NDVIA is:

[0051]

[0052] Among them, NDVIA i,j is the NDVI anomaly index of the jth month of the i-th year, i is the year sequence, j is the jth period of the growing season, the growing season is from April to September, and the NDVI ave,j is the multi-year NDVI mean of the jth period.

[0053] After calculating the monthly NDVIA, the annual NDVIA was then obtained by calculating the monthly average (considering only the growing season, i.e., April to September).

[0054] Execute step 130 to obtain MOD12Q1 products and extract land cover datasets for each year;

[0055] Execute step 140, use the optimal threshold extraction method to extract the optimal threshold for defining the drought index pixel by pixel, the optimal threshold extraction method is to calculate the correlation between the multi-year time series SMSDI data and the NDVIA data under different threshold conditions (the 46 thresholds corresponding to the 46 time series generated in step 110), and take the threshold corresponding to the maximum value of the correlation between the two as the optimal threshold for defining the drought index. The calculation formula of the correlation r is:

[0056]

[0057] Among them, SMSDI i SMSDI, NDVIA for the i-th year i is the NDVIA for the i-th year; is the mean value of SMSDI, are the mean values ​​of NDVIA; r max is the maximum correlation coefficient, and n is the total number of years.

[0058] The maximum correlation coefficient r max The calculation formula is

[0059] r max=max{r T0＝5th ,r T0＝6th ,r T0＝7th ,…,r T0＝50th}

[0060] Among them, r T0＝5th ,r T0＝6th ,r T0＝7th ,…,r T0＝50th is the correlation coefficient between SMSDI and NDVIA calculated when T0 takes the 5th, 6th, 7th, …, 50th percentile values ​​respectively, and max is the function of taking the maximum value in the set.

[0061] Execute step 150 to obtain the optimal threshold for defining the agricultural drought index under different underlying surfaces, including the following sub-steps: execute step 151, use the optimal threshold extraction method to determine the optimal threshold for each pixel; execute step 152, combine land cover data, and calculate the mean of the optimal threshold of all pixels under a specific underlying surface; execute step 153, use the mean to determine the optimal threshold for defining the drought index of the specific underlying surface.

[0062] Embodiment 2

[0063] The present invention aims to determine an index construction method that takes spatial heterogeneity into account. By this method, a drought index is constructed based on time series soil moisture and vegetation index data, so that the index has a significant agricultural drought indication effect.

[0064] In order to achieve the above object, the present invention provides the following scheme:

[0065] S1 first constructs the Generalized Soil-moisture-based Soil Dryness Index (GSMSDI) based on the time series soil moisture. The drought index is defined as follows:

[0066]

[0067] Where P sm is the percentile value of soil moisture in the time series data; T0 is the percentile threshold. When the soil moisture percentile value is lower than the threshold, GSMSDI < 0; when the soil moisture percentile value is higher than the threshold, GSMSDI > 0; otherwise, GSMSDI = 0. The GSMSDI index ranges from -1 to 1. The closer it is to 1, the wetter it is; the closer it is to -1, the drier it is. This generalized agricultural drought index has an uncertain drought threshold.

[0068] S2 defines the formula for the generalized agricultural drought index mentioned in S1, selects the value of T0 from 5 to 50, takes a value every step of 1, and finally determines the agricultural drought index (Soil-moisture-based Soil Dryness Index, SMSDI) corresponding to different thresholds. A total of 46 time series (T0 = 5, 6, 7, ..., 50th) are generated, such as Figure 2 As shown, 5 time series are selected, namely T0=10, 20, 30, 40, 50th. It can be seen that there are differences in the drought index time series obtained by different thresholds. Determining different thresholds has a greater impact on the judgment of "when to enter the drought state" (SMSDI<0)).

[0069] S3 obtains MODIS13A3 products and extracts the monthly Normalized Difference Vegetation Index (NDVI) data. Then, the NDVI anomaly (NDVIA) is calculated using the following formula using the NDVI data.

[0070]

[0071] Where i is the year sequence; j is the jth period (month) of the growing season (April to September); NDVIA i,j is the NDVI anomaly index of the jth period in the i-th year. ave,j The multi-year NDVI mean value for the jth period. The NDVI anomaly index obtained by formula (2) can eliminate the influence of vegetation growth during the growing season, so that the index can be compared over time. After calculating the monthly NDVIA, the annual NDVIA is obtained by calculating the monthly average value (considering only the growing season, i.e., April to September).

[0072] S4 obtains MOD12Q1 products and extracts land cover datasets for each year.

[0073] S5 Extract the optimal threshold value for defining drought index pixel by pixel: Calculate the correlation between SMSDI data and NDVIA data of multiple years, and take the threshold value corresponding to the maximum value of the correlation between the two as the optimal threshold value for defining the index. Figure 3 As shown in the figure, a pixel with grassland as the underlying surface is selected. The correlation between SMSDI and NDVIA shows a unimodal change characteristic with the change of threshold. When the threshold T0=25th, the correlation between drought index and vegetation growth condition is the strongest. The threshold at this time is defined as the optimal threshold. This method ensures that the drought index determined by the selected threshold has the greatest impact on vegetation.

[0074] S6 obtains the optimal threshold for defining the agricultural drought index under different underlying surfaces. Since the optimal threshold is determined at a regional scale, the heterogeneity of each pixel is relatively large. This method determines the optimal threshold pixel by pixel using the method mentioned in S5, and then calculates the average of the optimal thresholds of all pixels under a specific underlying surface type. Figure 4 The figure shows the statistical results of Inner Mongolia. The statistical results show that the optimal threshold of most pixels under grassland type is 10-19th, and the mean is 18th. The optimal threshold of drought index definition corresponding to grassland underlying surface is 18th. The mean is determined as the optimal threshold of drought index definition for a specific underlying surface, reflecting the spatial heterogeneity.

[0075] Verification and comparison of S7 results. Figure 5 As shown, three sites were selected and compared with the widely used SPEI index. It was found that SMSDI and SPEI had a more consistent trend of change and a better correlation. At the same time, the index has a strong correlation with the vegetation growth condition and is simple to calculate. In order to better understand the present invention, the above is described in detail in conjunction with the specific embodiments of the present invention, but it is not a limitation of the present invention. Any simple modification made to the above embodiments based on the technical essence of the present invention still belongs to the scope of the technical solution of the present invention. Each embodiment in this specification focuses on the differences from other embodiments, and the same or similar parts between the various embodiments can be referred to each other. For the system embodiment, since it basically corresponds to the method embodiment, the description is relatively simple, and the relevant parts can be referred to the partial description of the method embodiment.

Claims

1. A method for constructing an agricultural drought index considering spatial heterogeneity, comprising constructing a generalized agricultural drought index (GSMSDI), characterized in that: The following steps are also included: Step 1: Determine the agricultural drought index SMSDI corresponding to different thresholds, including constructing the generalized agricultural drought index GSMSDI based on the time series soil moisture. The generalized agricultural drought index GSMSDI is defined as Among them, P sm is the percentile value of soil moisture in the time series data, and T0 is the percentile threshold; For the generalized agricultural drought index GSMSDI, select T0 = 5, 6, 7, ..., 50th, take values ​​every step of 1, determine the agricultural drought index SMSDI corresponding to different thresholds, and generate a total of 46 time series; Step 2: Obtain MODIS13A3 products, extract monthly normalized difference vegetation index NDVI data, and use NDVI data to calculate NDVI anomaly index NDVIA. The formula of NDVI anomaly index NDVIA is: Among them, NDVIA i,j is the NDVI anomaly index of the jth month of the i-th year, i is the year sequence, j is the jth period of the growing season, the growing season is from April to September, and the NDVI ave,j is the multi-year NDVI mean of the jth period; After calculating the monthly NDVIA, the annual NDVIA is obtained by calculating the monthly average; Step 3: Obtain MOD12Q1 products and extract land cover datasets for each year; Step 4: Use the optimal threshold extraction method to extract the optimal threshold for defining the agricultural drought index pixel by pixel. The optimal threshold extraction method is to calculate the correlation between the multi-year time series SMSDI data and the NDVIA data under different threshold conditions, and the threshold corresponding to the maximum value of the correlation between the two is used as the optimal threshold for defining the agricultural drought index; Step 5: Obtaining the optimal threshold for defining the agricultural drought index under different underlying surfaces, including: using the optimal threshold extraction method to determine the optimal threshold pixel by pixel; combining land cover data to calculate the mean of the optimal threshold for all pixels under a specific underlying surface; and using the mean to determine the optimal threshold for defining the agricultural drought index for the specific underlying surface.

2. The method for constructing an agricultural drought index considering spatial heterogeneity according to claim 1, characterized in that: When the soil moisture percentile value is lower than T0, GSMSDI<0, when the soil moisture percentile value is higher than T0, GSMSDI>0, otherwise, GSMSDI=0.

3. The method for constructing an agricultural drought index considering spatial heterogeneity according to claim 2, characterized in that: The GSMSDI index ranges from -1 to 1. The closer it is to 1, the wetter it is, and the closer it is to -1, the drier it is.

Images