A method for fitting and reconstructing time series data function of vegetation index based on curve feature weighting

Abstract

The present invention discloses a method for fitting and reconstructing time series data of vegetation indices based on curve feature weighting. By utilizing the curve feature differences between noise points and vegetation signals in the time series curve of vegetation indices, high-quality time series data of vegetation indices are reconstructed through inputting remote sensing vegetation index time series data for a complete vegetation growth season, determining the gradual change points and sudden drop points in the vegetation index time series data, calculating the weights of each data point, and performing function weighted fitting. The characteristics of this method are that it does not rely on auxiliary data, has strong anti-noise ability and good universality, and can be used for reconstructing time series data of vegetation indices in regions with obvious seasonal changes of global vegetation.

Description

Technical Field

[0001] The invention relates to the field of vegetation remote sensing, and in particular to a remote sensing vegetation index time series data reconstruction method. Background Art

[0002] Vegetation index time series data obtained by satellite sensors can better reflect the growth status and seasonal changes of vegetation, and are widely used in land surface vegetation parameter estimation and land cover classification. However, the quality of vegetation index time series data is easily affected by clouds and atmospheric conditions, which will cause errors in subsequent applications. Although most vegetation index time series data products are generated using the maximum synthesis method and a series of cloud detection algorithms, they still contain a lot of residual noise.

[0003] At present, a series of vegetation index time series data reconstruction methods have been developed, which can be roughly divided into four categories: (1) time domain filtering methods, such as Savitzky-Golay filtering method, mean iteration filtering method and optimal index slope extraction method; (2) time domain function fitting methods, such as dual logistic function fitting method, dual Gaussian function fitting method and polynomial function fitting method; (3) frequency domain filtering methods, such as time series harmonic analysis method; (4) time-space fusion methods, such as time-space fusion Savitzky-Golay filtering method. The key parameter settings of time domain filtering method and frequency domain filtering method usually need to be subjectively set according to the researcher's experience, which will lead to new uncertainties; time domain function fitting method is easily affected by noise points in vegetation index time series data; time-space fusion method is easy to fail when applied to images with large cloud cover.

[0004] Existing studies often use additional auxiliary data to reduce the impact of noise on the vegetation index time series data reconstructed by the function fitting method. Auxiliary data are usually used in two ways: (1) based on the cloud markers in the auxiliary data, the vegetation index points marked as clouds are removed and replaced with the linear interpolation results of adjacent vegetation index points not affected by clouds; (2) based on the signal-to-noise ratio markers in the auxiliary data, the fitting weights of each point in the vegetation index time series data are set. Although these two methods can effectively improve the quality of vegetation index time series data reconstruction by the function fitting method, the function fitting method based on auxiliary data has great limitations due to the limited accuracy of the auxiliary data itself, the poor universality of the auxiliary data processing rules, and the lack of effective auxiliary data in some data sources.

[0005] In fact, the time-series data of vegetation indices itself contains the quality information of each data point. Vegetation index points with high signal-to-noise ratio reflect the gradual process of vegetation growth and senescence, while vegetation index points with a sudden drop are usually contaminated by clouds and adverse atmospheric conditions. We can use this quality information to weight each data point of the time-series data of vegetation indices, thereby effectively avoiding the dependence on auxiliary data in the function fitting method for reconstructing the time-series data of vegetation indices. Summary of the Invention

[0006] In view of some deficiencies of the traditional function fitting reconstruction method based on auxiliary data, the present invention proposes a function fitting reconstruction method for time-series data of vegetation indices based on curve feature weighting, aiming to weight each data point using the curve features of the time-series data of vegetation indices and reconstruct high-quality time-series data of vegetation indices.

[0007] To achieve this purpose, the present invention adopts the following technical solutions:

[0008] A function fitting reconstruction method for time-series data of vegetation indices based on curve feature weighting, comprising the following steps:

[0009] A. Input the time-series data of remote sensing vegetation indices for a complete vegetation growing season.

[0010] B. Determine the gradual change points and sudden drop points in the time-series data of vegetation indices. All data points that show a monotonically increasing trend from both ends to the peak point in the time-series data of vegetation indices are determined as gradual change points, and the remaining data points are determined as sudden drop points. The specific steps for determining the gradual change points and sudden drop points in the time-series data of vegetation indices are as follows:

[0011] (1) Use the time position of the peak point of the vegetation index to divide the time-series curve of the vegetation index into two parts, the left and the right, and then determine the gradual change points and sudden drop points respectively.

[0012] (2) In the rising stage of the left curve, judge sequentially from the first point to the peak point. If the vegetation index value of a certain point is greater than or equal to the vegetation index values of all points before this point, then this point is a gradual change point; otherwise, it is a sudden drop point.

[0013] (3) In the falling stage of the right curve, judge sequentially from the last point to the peak point. If the vegetation index value of a certain point is greater than or equal to the vegetation index values of all points after this point, then this point is a gradual change point; otherwise, it is a sudden drop point.

[0014] C. Calculate the weights. The weights of the gradual change points are uniformly set to 1, and the weights of the sudden drop points are comprehensively determined according to their drop amplitude and proximity to the peak. The specific steps for weight calculation are as follows:

[0015] (1) The weights of the gradual change points in the time-series data of vegetation indices are uniformly set to 1, and its mathematical model is:

[0016] W GCP = 1

[0017] (2) To determine the weight of the sudden drop points in the time series data of the vegetation index, first linearly stretch the time series data of the vegetation index to 0 - 10, and its mathematical model is:

[0018]

[0019] (3) The drop amplitude of the sudden drop point is the drop amplitude of its stretched VI Ls value, and its mathematical model is:

[0020]

[0021] (4) The peak proximity of the sudden drop point is the proximity degree in time between the sudden drop point and the peak point of the vegetation index, and its mathematical model is:

[0022]

[0023] (5) The weight of the sudden drop point is comprehensively determined according to its drop amplitude and peak proximity, and its mathematical model is:

[0024]

[0025] In the above mathematical model, W GCP is the weight of the gradual change point; VI is the original vegetation index value, and VI Ls represents the result of linearly stretching the time series data of the vegetation index to 0 - 10, and VI max and VI min are the maximum and minimum values of the original time series data of the vegetation index respectively; Δh is the drop amplitude of the sudden drop point, and VI LS_SDP , VI Ls_GCP1 , VI LS_GCP2 are the VI LS values of the sudden drop point and its adjacent front and rear gradual change points respectively, and DOY SDP , DOY GCP1 , DOY GCP2 are the day of year of the sudden drop point and its adjacent front and rear gradual change points respectively; P is the peak proximity of the sudden drop point, and DOY m is the day of year of the peak point of the vegetation index, and DOY min and DOY max are the maximum and minimum values of the day of year in the time series data of the vegetation index respectively; W SDP is the weight of the sudden drop point; if the day of year in the formula crosses to the next year, the day of year of all data points in the next year needs to be added with 365.

[0026] D. Perform function weighted fitting. Based on the weights of each data point that have been determined, use an S-shaped function or a polynomial function to perform function weighted fitting and reconstruction on the time series data of the vegetation index. Among them, the S-shaped functions include double logistic functions, double Gaussian functions, etc.

[0027] The present invention has the following characteristics:

[0028] (1) It does not rely on auxiliary data, and weights each data point by using the curve characteristics of the time series data of the vegetation index;

[0029] (2) It is simple and easy to implement, and the reconstruction result is reliable;

[0030] (3) It has a wide range of applicability and can be applied to the reconstruction of the time series data of the vegetation index in regions with obvious seasonal changes of global vegetation. Description of the Drawings

[0031] Figure 1 It is the reconstruction result (upper) of the time series data of the vegetation index and the fitting weights (lower) of each point. Detailed Embodiment

[0032] The following further describes the technical implementation scheme of the present invention with reference to the drawings.

[0033] Based on the MOD09A1 data (time resolution of 8 days and spatial resolution of 500m), calculate the NDVI (Normalized Difference Vegetation Index) time series data of a grassland pixel (34.33°N, 101.72°E) on the Qinghai-Tibet Plateau in 2018, and then use the method in the present invention to reconstruct it.

[0034] A. Input the time series data of the remote sensing vegetation index for a complete vegetation growth season.

[0035] In this case, the vegetation type of the pixel is grassland, its growth cycle is one year, and the NDVI time series data in 2018 is for a complete vegetation growth season, with a total of 46 data points in a year.

[0036] B. Determine the gradual change points and sudden drop points in the time series data of the vegetation index. Determine all the data points that maintain a monotonically increasing trend from both ends to the peak point in the time series data of the vegetation index as the gradual change points, and the remaining data points as the sudden drop points. The specific steps to determine the gradual change points and sudden drop points in the NDVI time series data are as follows:

[0037] (1) Use the time position where the NDVI peak point is located to divide the NDVI time series curve into two parts on the left and right, and then determine the gradual change points and sudden drop points respectively.

[0038] (2) In the rising stage of the left curve, judge successively from the first point to the peak point. If the NDVI value of a certain point is greater than or equal to the NDVI values of all the points before this point, then this point is a gradual change point; otherwise, it is a sudden drop point.

[0039] (3) In the falling stage of the right curve, judge successively from the last point to the peak point. If the NDVI value of a certain point is greater than or equal to the NDVI values of all the points after this point, then this point is a gradual change point; otherwise, it is a sudden drop point.

[0040] C. Calculate the weights. The weights of the gradual change points are uniformly set to 1, and the weights of the sudden drop points are comprehensively determined according to their drop amplitudes and peak proximities. The specific steps for weight calculation are as follows:

[0041] (1) The weights of the gradual change points in the NDVI time series data are uniformly set to 1, and its mathematical model is:

[0042] W GCP = 1

[0043] (2) To determine the weights of the sudden drop points in the NDVI time series data, first linearly stretch the NDVI time series data to 0 - 10, and its mathematical model is:

[0044]

[0045] (3) The drop amplitude of the sudden drop point is the drop amplitude of its stretched NDVI LS value, and its mathematical model is:

[0046]

[0047] (4) The peak proximity of the sudden drop point is the proximity degree of the sudden drop point and the NDVI peak point in time, and its mathematical model is:

[0048]

[0049] (5) The weights of the sudden drop points are comprehensively determined according to their drop amplitudes and peak proximities, and its mathematical model is:

[0050]

[0051] In the above mathematical models, W GCP is the weight of the gradual change point; NDVI is the NDVI value, NDVI LS represents the result of linearly stretching the NDVI time series data to 0 - 10, NDVI max and NDVI min are respectively the maximum value and the minimum value of the NDVI time series data; Δh is the drop amplitude of the sudden drop point, NDVI LS_SDP , NDVI LS_GCP1, NDVI LS_GCP2 are the NDVI values of the sudden drop point and its adjacent gradual change points before and after respectively LS value, DOY SDP , DOY GCP1 , DOY GCP2 are the annual day sequences of the sudden drop point and its adjacent gradual change points before and after respectively; P is the peak proximity of the sudden drop point, and DOY m is the annual day sequence of the NDVI peak point, DOY min and DOY max are the maximum and minimum values of the annual day sequence in the NDVI time series data respectively; W SDP is the weight of the sudden drop point; if the annual day sequence in the formula spans to the next year, 365 needs to be added to the annual day sequences of all data points in the next year.

[0052] D. Perform function weighted fitting. Based on the weights of each determined data point, use the S-shaped function or polynomial function to perform function weighted fitting and reconstruction on the vegetation index time series data. Among them, the S-shaped functions include double logistic functions, double Gaussian functions, etc.

[0053] In this case, based on the weights of each determined data point, select the double logistic function to perform function weighted fitting and reconstruction on the NDVI time series data in 2018.

[0054] As mentioned above, it is only the preferred specific implementation manner of the present invention, but the protection scope of the present invention is not limited thereto. Any changes or substitutions that can be easily thought of by those familiar with the technology within the technical scope disclosed by the present invention should be covered within the protection scope of the present invention. Therefore, the protection scope of the present invention should be subject to the protection scope of the claims.

Claims

1. A method for function fitting and reconstruction of time series data of vegetation indices based on curve feature weighting, characterized in that, it includes the following steps: A. Input the time series data of remote sensing vegetation indices for a complete vegetation growth season; B. Determine the gradual change points and sudden drop points in the time series data of vegetation indices. All data points that maintain a monotonically increasing trend from both ends to the peak point in the time series data of vegetation indices are determined as gradual change points, and the remaining data points are determined as sudden drop points; C. Calculate the weights. The weights of the gradual change points are uniformly set to 1, and the weights of the sudden drop points are comprehensively determined according to their drop amplitude and peak proximity; D. Perform function weighted fitting. Based on the weights of each data point that have been determined, use an S-shaped function or a polynomial function to perform function weighted fitting and reconstruction on the time series data of vegetation indices; The method for determining the weights of the sudden drop points in step C is the following 3 mathematical formulas: Where Δh is the drop amplitude at the sudden drop point, VI LS represents the result of linearly stretching the time series data of the vegetation index to 0-10, VI LS_SDP , VI LS_GCP1 , VI LS_GCP2 are the VI values of the sudden drop point and its adjacent gradual change points before and after respectively, DOY LS , DOY SDP , DOY GCP1 , DOY GCP2 are the annual day sequences of the sudden drop point and its adjacent gradual change points before and after respectively, P is the peak proximity of the sudden drop point, DOY m is the annual day sequence of the maximum value point of the vegetation index, DOY min and DOY max are the maximum and minimum values of the annual day sequence in the time series data of the vegetation index respectively, W SDP is the weight of the sudden drop point; if the annual day sequence in the formula spans to the next year, 365 needs to be added to the annual day sequences of all data points in the next year.

Images