A remote sensing ecological index monitoring method and device suitable for long time series

Abstract

The application relates to a remote sensing ecological index monitoring method and device suitable for long time sequences. The method and device perform data synthesis on acquired long time sequence Landsat remote sensing data; calculate time sequence remote sensing indexes of the Landsat remote sensing data; detect invariable regions of the time sequence remote sensing indexes of the Landsat remote sensing data; normalize the time sequence remote sensing indexes in the probability space of the invariable regions; uniformly perform principal component analysis on normalized indexes of multiple time phases to acquire unified principal component load coefficients. The application is based on medium spatial resolution Landsat remote sensing images, utilizes time sequence data synthesis, time sequence invariable region detection, time sequence principal component analysis and other methods to design a remote sensing ecological index suitable for long time sequence dramatic change regions, and finally accurately reflects the distribution and evolution conditions of ecological environment conditions in a large area range and a long time sequence.

Description

Technical Field

[0001] This invention relates to the field of ecological remote sensing monitoring technology, and more specifically, to a method and apparatus for monitoring remote sensing ecological indices suitable for long-term time series. Background Technology

[0002] A healthy ecological environment is a crucial prerequisite for human survival and development. Human activities have a significant impact on the ecological environment, thus necessitating its monitoring and assessment. Currently, remote sensing technology, with its advantages of rapidly and accurately acquiring large-scale information, is widely used in various fields, such as ecological environment quality assessment.

[0003] However, methods using single remote sensing indices to assess the environment can only monitor one type of environmental characteristic, such as using vegetation indices to monitor forest cover, surface temperature to monitor the urban heat island effect, and building indices to monitor urban expansion rates. These methods cannot provide a comprehensive evaluation of complex ecological environments. Multi-indice fusion remote sensing ecological indices are an effective way to conduct comprehensive evaluations of complex ecological environments; however, they currently face three main challenges:

[0004] (1) In cloudy and rainy areas, taking South my country as an example, it is difficult to guarantee the acquisition of continuous and effective optical images at the same time when conducting long-term ecological monitoring.

[0005] (2) In areas with drastic changes in the ecological environment, taking the Guangdong-Hong Kong-Macao Greater Bay Area of ​​my country as an example, due to the difference in data time phases, the monitoring results of a single time phase mainly reflect the relative status of the region in a single time phase. However, when conducting long-term series monitoring, the monitoring results of multiple time phases often cannot accurately reflect the true dynamics of the time series.

[0006] (3) In large-area ecological monitoring, previous studies were mainly based on low-resolution MODIS images. When applied to urban areas, they could not describe the distribution and changes of ecological conditions in detail.

[0007] Therefore, these factors reduce the stability, accuracy, and wide applicability of remote sensing monitoring of the ecological environment.

[0008] Xu Hanqiu proposed a Remote Sensing Ecology Index (RSEI) that couples multiple indicators. This index couples four indicators: greenness, humidity, heat, and aridity. Principal component analysis is used to determine the weights, avoiding the uncertainty of manually determining the weights, and it has been widely used. Furthermore, the RSEI remote sensing ecological index is constructed based on the Stress-State-Response (PSR) framework. This framework links the environmental stresses, the current state of the ecological environment, and the response to changes in the ecological environment, evaluating changes in the ecological environment under the influence of human factors and other influences.

[0009] The Remote Sensing Ecological Index (RSEI) integrates greenness, humidity, dryness, and heat. Dryness represents the current environmental stress, greenness represents the current state of the ecological environment, and humidity and heat describe the response of the local ecological environment to changes in climate conditions. Greenness is represented by the Normalized Difference Vegetation Index (NDVI), humidity by the humidity component in the tasseling transformation, heat by land surface temperature (LST), and dryness by the average of the Integral Building Index (IBI) and the Bare Soil Index (SI).

[0010] Since the dimensions of the aforementioned humidity, greenness, heat, and dryness indices are inconsistent, to facilitate calculation and ensure consistent weights for the four indices, they must first be normalized so that their ranges are all between [0, 1]. After normalization, the four normalization factors are combined into one image. Principal component analysis (PCA) is then performed using the tool Forward PCA Rotation New Statistics and Rotate to extract the first principal component and calculate the initial RSEI value, denoted as RSEI0.

[0011] RSEI0=1-PC1[f(NDVI, Wet, LST, NDBSI)]

[0012] The RSEI values ​​vary across different periods and regions. To facilitate comparative studies, RSEI0 is normalized again:

[0013] RSEI = (RSEI0 - RSEI) 0_min ) / (RSEI 0_max -RSEI 0_min )

[0014] The final RSEI is obtained from the above formula. Its value ranges between [0, 1]. The higher the RSEI value, the better the ecological quality, and vice versa.

[0015] However, the aforementioned existing technologies have three main problems:

[0016] (1) In the RSEI inversion of a single time phase, the water body has a high humidity index but a low greenness index, resulting in a low RSEI value for the water body in the final RSEI, which means that the ecological condition is poor, which is quite different from people's general understanding.

[0017] (2) In the RSEI inversion of a single time phase, the four indicators need to be normalized first. However, the normalization results based on the maximum and minimum values ​​are often not good, especially for the thermal index. Using the extreme high temperature of some artificial surfaces (such as landfills) as the maximum value for normalization will result in the normalized thermal index of the entire study area being too low, which will affect the final RSEI results.

[0018] (3) In the multi-temporal RSEI inversion, the relationship between the temporal phases is not considered. The RSEI inversion of each temporal phase is based on a single temporal phase. The final RSEI result reflects the relative situation in a single temporal phase area. When the ecological situation in the study area changes drastically between multiple temporal phases, the temporal change characteristics of the multi-temporal RSEI result of the same pixel cannot be consistent with the actual ecological change situation. Summary of the Invention

[0019] This invention provides a method and apparatus for monitoring remote sensing ecological indices over long time series, which at least solves the technical problem of low accuracy in existing remote sensing ecological index monitoring results.

[0020] According to an embodiment of the present invention, a method for monitoring remote sensing ecological indices suitable for long-term time series is provided, comprising the following steps:

[0021] Data synthesis was performed on the acquired long-term Landsat remote sensing data;

[0022] Calculate the time-series remote sensing index of the Landsat remote sensing data after data synthesis;

[0023] Mann-Kendall trend test was performed on the time series remote sensing indices of Landsat remote sensing data to detect invariant regions.

[0024] Normalization of time series remote sensing indices is performed using the probability space of invariant regions;

[0025] By unifying the normalized indices of multiple time phases and performing principal component analysis, a unified principal component loading coefficient is obtained.

[0026] Furthermore, the data synthesis of the acquired long-term Landsat remote sensing data includes:

[0027] Input Landsat time-series image data of the same region, perform cloud cover detection on each image, mask the cloud cover pixels, and synthesize the data by selecting the nearest N years of peak vegetation growth period. The calculation formula is as follows:

[0028]

[0029] Where y represents the year of synthesis, and y-1≤yi ≤y+1, d1≤d j ≤d2, median is the median operator.

[0030] Furthermore, after synthesizing the acquired long-term Landsat remote sensing data, the method also includes:

[0031] For each synthesized Landsat image, water bodies are detected, and the pixels covered by water bodies are masked.

[0032] Furthermore, the time-series remote sensing indices of the synthesized Landsat remote sensing data are calculated as follows:

[0033] For each composite Landsat image, the remote sensing indices NDVI, NDWI, IBI, SI, and temperature LST were calculated using the following formulas:

[0034] NDVI = (B nir -B red ) / (B nir +B red )

[0035] NDWI = (B green -B nir ) / (B green +B nir )

[0036]

[0037]

[0038]

[0039] in L6 = gain × DN + bias;

[0040] L6 is the radiometrically calibrated thermal infrared reflectivity. For the TM sensor, K1 = 607.76 W·m⁻²·sr⁻¹·μm⁻¹, K2 = 1260.56 K; for the TIRS sensor in band 10, K1 = 774.89 W·m⁻²·sr⁻¹·μm⁻¹, K2 = 1321.08 K; λ is the center wavelength of the thermal infrared band. TM =11.435μm, λ b10 =10.900μm; ρ =1.438×10⁻²m·K; ε is the surface emissivity;

[0041] Based on the tassel transformation, the humidity index Wet is established as follows:

[0042] Wet TM=0.0315×B blue +0.2021×B green +0.3012×Bρ red +0.1594×B nir -0.6806×B swir1 -0.6109B swir2

[0043] Wet OLI =0.1511×B blue +0.1973×B green +0.3283×B red +0.3407×B nir -0.7117×B swir1 -0.4559B swir2

[0044] Wet TM Wet OLI These correspond to the Landsat5TM sensor and the Landsat8OLI sensor, respectively.

[0045] Based on IBI and SI, the dryness index NDBSI is established, and the formula is as follows:

[0046]

[0047] NDVI was used as the greenness index, and LST was used as the heat index.

[0048] Furthermore, the Mann-Kendall trend test was performed on the time-series remote sensing indices of Landsat remote sensing data to detect invariant regions, including:

[0049] Mann-Kendall trend tests were performed on time series data for NDVI, NDWI, IBI, and SI respectively, and invariant layers were created. Taking NDVI as an example, layers were created. The calculation formula is as follows:

[0050] UC = Mann-Kendall{NDVI i}, i = the time series being monitored.

[0051]

[0052] Compare with the above formula and establish respectively And create an invariant layer. UC The formula is as follows:

[0053]

[0054] Furthermore, normalization of time-series remote sensing indices using the probability space of invariant regions includes:

[0055] Based on invariant layer UC Time series data for NDVI (Greenness Dispersion), Wet (Humidity), NDBSI (Dryness Dispersion), and LST (Heat Dispersion) were masked. For each layer of the masked time phase, probability space statistics ranging from 2% to 98% were performed, with 2% as the minimum and 98% as the maximum. The original unmasked greenness, humidity, dryness, and heat layers were then normalized. Taking 2020 NDVI as an example, the normalized greenness layer is as follows: The calculation formula is as follows:

[0056]

[0057] The dryness and heat indices are inversely normalized. Taking heat LST as an example, the normalization is as follows:

[0058]

[0059] Furthermore, the normalized indices from multiple time phases are uniformly analyzed using principal component analysis to obtain unified principal component loading coefficients, including:

[0060] Several sample points were randomly deployed within the study area, with the number of sample points determined based on the study area's extent and data resolution. Based on these sample points, four normalized indicators from multiple time phases were extracted annually. These four normalized indicators from the time series were then combined into four new indicator columns. and Based on the above index list, PCA principal component analysis was performed to obtain the loading coefficients of the first principal component, labeled as a, b, c, and d, respectively, corresponding to... and Four indicator columns.

[0061] Furthermore, after uniformly performing principal component analysis on the normalized indices of multiple time phases to obtain unified principal component loading coefficients, the method also includes:

[0062] For each phase of the observed time series, the initial RSEI0 value is calculated using the following formula:

[0063]

[0064] Where y represents the monitoring time phase;

[0065] Acquired for all observation phases Statistical maximum value RSEI across multiple time phases max and minimum value RSEI minThen, normalize the result to obtain the final RSEI, as shown in the following formula:

[0066]

[0067] Furthermore, after uniformly performing principal component analysis on the normalized indices of multiple time phases to obtain unified principal component loading coefficients, the method also includes:

[0068] The regions of the yearly cloud masks, after masking the water-covered pixels, were applied to all observation phases. Statistical maximum value RSEI across multiple time phases max and minimum value RSEI min Then, normalize it to obtain the final RSEI. This is part of the RSEI generation process. y And assign the cloud mask area an empirical value.

[0069] According to another embodiment of the present invention, a remote sensing ecological index monitoring device suitable for long-term time series is provided, comprising:

[0070] The data synthesis unit is used to synthesize the acquired long-term Landsat remote sensing data.

[0071] The remote sensing index calculation unit is used to calculate the time series remote sensing index of the Landsat remote sensing data after data synthesis.

[0072] The trend test unit is used to perform the Mann-Kendall trend test on the time series remote sensing indices of Landsat remote sensing data to detect invariant regions.

[0073] The normalization unit is used to normalize the time series remote sensing index using the probability space of the invariant region;

[0074] The principal component analysis unit is used to perform principal component analysis on normalized indices from multiple time phases to obtain unified principal component loading coefficients.

[0075] A storage medium storing program files capable of implementing any of the above-mentioned remote sensing ecological index monitoring methods applicable to long-term time series.

[0076] A processor for running a program, wherein the program executes any of the above-mentioned remote sensing ecological index monitoring methods applicable to long-term series.

[0077] The remote sensing ecological index monitoring method and device applicable to long-term series of the present invention are based on medium spatial resolution Landsat remote sensing images. They utilize methods such as time series data synthesis, time series invariant region detection, and time series principal component analysis to design a remote sensing ecological index applicable to regions with drastic changes over a long period of time. Ultimately, it accurately reflects the distribution and evolution of the ecological environment over a large area and a long time series. After extensive experimental verification, the remote sensing ecological index monitoring results of the present invention have high accuracy, the operation process is simple and feasible, and it has considerable practical value. Attached Figure Description

[0078] The accompanying drawings, which are included to provide a further understanding of the invention and form part of this application, illustrate exemplary embodiments of the invention and, together with their description, serve to explain the invention and do not constitute an undue limitation thereof. In the drawings:

[0079] Figure 1 This is a flowchart of the remote sensing ecological index monitoring method applicable to long-term time series of the present invention;

[0080] Figure 2 This is a schematic diagram of the remote sensing ecological index monitoring results of the Guangdong-Hong Kong-Macao Greater Bay Area from 1990 to 2020 based on the technical solution of this invention.

[0081] Figure 3 This is a block diagram of the remote sensing ecological index monitoring device applicable to long-term time series according to the present invention. Detailed Implementation

[0082] To enable those skilled in the art to better understand the present invention, the technical solutions of the present invention will be clearly and completely described below with reference to the accompanying drawings of the embodiments of the present invention. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort should fall within the scope of protection of the present invention.

[0083] It should be noted that the terms "first," "second," etc., in the specification, claims, and accompanying drawings of this invention are used to distinguish similar objects and are not necessarily used to describe a specific order or sequence. It should be understood that such data can be interchanged where appropriate so that the embodiments of the invention described herein can be implemented in orders other than those illustrated or described herein. Furthermore, the terms "comprising" and "having," and any variations thereof, are intended to cover a non-exclusive inclusion; for example, a process, method, system, product, or apparatus that comprises a series of steps or units is not necessarily limited to those steps or units explicitly listed, but may include other steps or units not explicitly listed or inherent to such processes, methods, products, or apparatus.

[0084] Example 1

[0085] According to an embodiment of the present invention, a remote sensing ecological index monitoring method suitable for long-term time series is provided, see [link to relevant documentation]. Figure 1 This includes the following steps:

[0086] S100: Perform data synthesis on the acquired long-term Landsat remote sensing data;

[0087] S200: Calculate the time series remote sensing index of the Landsat remote sensing data after data synthesis;

[0088] S300: Use the Mann-Kendall trend test on the time-series remote sensing indices of Landsat remote sensing data to detect invariant regions;

[0089] S400: Normalization of time series remote sensing indices using the probability space of invariant regions;

[0090] S500: Performs principal component analysis on normalized indices from multiple time phases to obtain unified principal component loading coefficients.

[0091] The remote sensing ecological index monitoring method applicable to long-term time series in this embodiment of the invention is based on medium spatial resolution Landsat remote sensing imagery. It uses methods such as time series data synthesis, time series invariant region detection, and time series principal component analysis to design a remote sensing ecological index applicable to regions with drastic changes over a long period of time. Ultimately, it accurately reflects the distribution and evolution of the ecological environment over a large area and a long time series. After extensive experimental verification, the remote sensing ecological index monitoring results of this invention have high accuracy, the operation process is simple and feasible, and it has considerable practical value.

[0092] The data synthesis of the acquired long-term Landsat remote sensing data includes:

[0093] Input Landsat time-series image data of the same region, perform cloud cover detection on each image, mask the cloud cover pixels, and synthesize the data by selecting the nearest N years of peak vegetation growth period. The calculation formula is as follows:

[0094]

[0095] Where y represents the year of synthesis, and y-1≤y i ≤y+1, d1≤d j ≤d2, median is the median operator.

[0096] The method further includes, after data synthesis of the acquired long-term Landsat remote sensing data, the following steps:

[0097] For each synthesized Landsat image, water bodies are detected, and the pixels covered by water bodies are masked.

[0098] The time-series remote sensing indices of the synthesized Landsat remote sensing data include:

[0099] For each composite Landsat image, the remote sensing indices NDVI, NDWI, IBI, SI, and temperature LST were calculated using the following formulas:

[0100] NDVI = (B nir -B red ) / (B nir +B red )

[0101] NDWI = (B green -B nir ) / (B green +B nir )

[0102]

[0103]

[0104]

[0105] in L6 = gain × DN + bias;

[0106] L6 is the radiometrically calibrated thermal infrared reflectivity. For the TM sensor, K1 = 607.76 W·m⁻²·sr⁻¹·μm⁻¹, K2 = 1260.56 K; for the TIRS sensor in band 10, K1 = 774.89 W·m⁻²·sr⁻¹·μm⁻¹, K2 = 1321.08 K; λ is the center wavelength of the thermal infrared band. TM =11.435μm, λ b10 =10.900μm; ρ =1.438×10⁻²m·K; ε is the surface emissivity;

[0107] Based on the tassel transformation, the humidity index Wet is established as follows:

[0108] Wet TM =0.0315×B blue +0.2021×B green +0.3012×Bρ red +0.1594×Bnir -0.6806×B swir1 -0.6109B swir2

[0109] Wet OLI =0.1511×B blue +0.1973×B green +0.3283×B red +0.3407×B nir -0.7117×B swir1 -0.4559B swir2

[0110] Wet TM Wet OLI These correspond to the Landsat5TM sensor and the Landsat8OLI sensor, respectively.

[0111] Based on IBI and SI, the dryness index NDBSI is established, and the formula is as follows:

[0112]

[0113] NDVI was used as the greenness index, and LST was used as the heat index.

[0114] Among these methods, the Mann-Kendall trend test is performed on the time-series remote sensing indices of Landsat remote sensing data to detect invariant regions, including:

[0115] Mann-Kendall trend tests were performed on time series data for NDVI, NDWI, IBI, and SI respectively, and invariant layers were created. Taking NDVI as an example, layers were created. The calculation formula is as follows:

[0116] UC = Mann-Kendall{NDVI i}, i = the time series being monitored.

[0117]

[0118] Compare with the above formula and establish respectively And create an invariant layer. UC The formula is as follows:

[0119]

[0120] The normalization of time-series remote sensing indices using the probability space of invariant regions includes:

[0121] Based on invariant layer UCTime series data for NDVI (Greenness Dispersion), Wet (Humidity), NDBSI (Dryness Dispersion), and LST (Heat Dispersion) were masked. For each layer of the masked time phase, probability space statistics ranging from 2% to 98% were performed, with 2% as the minimum and 98% as the maximum. The original unmasked greenness, humidity, dryness, and heat layers were then normalized. Taking 2020 NDVI as an example, the normalized greenness layer is as follows: The calculation formula is as follows:

[0122]

[0123] The dryness and heat indices are inversely normalized. Taking heat LST as an example, the normalization is as follows:

[0124]

[0125] This involves unifying the normalized indices from multiple time phases through principal component analysis to obtain unified principal component loading coefficients, including:

[0126] Several sample points were randomly deployed within the study area, with the number of sample points determined based on the study area's extent and data resolution. Based on these sample points, four normalized indicators from multiple time phases were extracted annually. These four normalized indicators from the time series were then combined into four new indicator columns. and Based on the above index list, PCA principal component analysis was performed to obtain the loading coefficients of the first principal component, labeled as a, b, c, and d, respectively, corresponding to... and Four indicator columns.

[0127] The method further includes, after uniformly performing principal component analysis on the normalized indices of multiple time phases to obtain uniform principal component loading coefficients, the following:

[0128] For each phase of the observed time series, the initial RSEI0 value is calculated using the following formula:

[0129]

[0130] Where y represents the monitoring time phase;

[0131] Acquired for all observation phases Statistical maximum value RSEI across multiple time phases max and minimum value RSEI min Then, normalize the result to obtain the final RSEI, as shown in the following formula:

[0132]

[0133] The method further includes, after uniformly performing principal component analysis on the normalized indices of multiple time phases to obtain uniform principal component loading coefficients, the following:

[0134] The regions of the yearly cloud masks, after masking the water-covered pixels, were applied to all observation phases. Statistical maximum value RSEI across multiple time phases max and minimum value RSEI min Then, normalize it to obtain the final RSEI. This is part of the RSEI generation process. y And assign the cloud mask area an empirical value.

[0135] The following detailed description, using specific embodiments, illustrates the remote sensing ecological index monitoring method of the present invention applicable to long-term time series:

[0136] To address the shortcomings of existing technologies, this invention tackles issues in current remote sensing ecological index monitoring of ecological conditions, such as single-temporal normalization and multi-temporal comparison. It utilizes long-term Landsat remote sensing data as the data source to establish a remote sensing ecological index monitoring method suitable for long-term, drastically changing regions. This invention uses regional vegetation growth windows for data synthesis and employs the Mann-Kendall trend test for time-series remote sensing indices (NDVI, NDWI, IBI, SI) to detect invariant areas. NDVI, humidity, IBI, SI, and LST are normalized using a probability space of 2%–98% for invariant areas to avoid distortion caused by extreme maximum or minimum values. In the principal component analysis step, the normalized indices from multiple temporal phases are uniformly analyzed to obtain unified principal component loading coefficients. For water bodies, water body detection is used to mask the water body, and in the final RSEI results, gaps in the water body data are filled by direct assignment.

[0137] This invention is applicable to the overall evaluation of regional ecological conditions in cloudy and rainy areas, covering a wide range and over a long period of time. For areas where the ecological conditions have changed drastically, it can accurately express both the relative characteristics of the regional ecological conditions in a single time phase and the temporal changes of the ecological conditions in multiple time phases. The comprehensive mapping results of multiple time phases can help analyze the spatiotemporal evolution of the regional ecological conditions over a long period of time.

[0138] Specifically, this invention, based on medium spatial resolution Landsat remote sensing imagery, utilizes methods such as time-series data synthesis, time-series invariant region detection, and time-series principal component analysis to design a remote sensing ecological index suitable for regions experiencing drastic long-term changes. Ultimately, it accurately reflects the distribution and evolution of the ecological environment over a large area and over a long time period. Extensive experimental verification has demonstrated that the remote sensing ecological index of this invention achieves high accuracy, has a simple and feasible operation process, and possesses considerable practical value. The detailed technical solution of this invention is described below:

[0139] A. Input time-series Landsat image data of the same region, perform cloud cover detection on each image, and mask the cloud cover pixels. Select the peak vegetation growth period of the past 3 years for data synthesis. Taking South China as an example, the peak growth period is selected as September-October. The calculation formula is as follows:

[0140]

[0141] Where y represents the year of synthesis, and y-1≤y i ≤y+1, d1≤d j ≤d2, where d1 represents September 1st, d2 represents October 31st, and median is the median operator.

[0142] B. Perform water body detection on each synthesized Landsat image, mask the water-covered pixels, and save the water mask file for each image.

[0143] C. For each composite Landsat image, calculate the remote sensing indices: NDVI, NDWI, IBI, SI, and temperature LST. The calculation formulas are as follows:

[0144] NDVI = (B nir -B red ) / (B nir +B red )

[0145] NDWI = (B green -B nir ) / (B green +B nir )

[0146]

[0147]

[0148]

[0149] in L6 = gain × DN + bias;

[0150] L6 is the radiometrically calibrated thermal infrared reflectivity. For the TM sensor, K1 = 607.76 W·m⁻²·sr⁻¹·μm⁻¹, K2 = 1260.56 K; for the TIRS sensor in band 10, K1 = 774.89 W·m⁻²·sr⁻¹·μm⁻¹, K2 = 1321.08 K; λ is the center wavelength of the thermal infrared band. TM =11.435μm, λ b10 =10.900μm; ρ =1.438×10-2m·K; ε is the surface emissivity.

[0151] Based on the tassel transformation, the humidity index Wet is established as follows:

[0152] Wet TM =0.0315×B blue +0.2021×B green +0.3012×Bρ red +0.1594×B nir -0.6806×B swir1 -0.6109B swir2

[0153] Wet OLI =0.1511×B blue +0.1973×B green +0.3283×B red +0.3407×B nir -0.7117×B swir1 -0.4559B swir2

[0154] Wet TM Wet OLI These correspond to the Landsat5TM sensor and the Landsat8OLI sensor, respectively.

[0155] Based on IBI and SI, the dryness index NDBSI is established, and the formula is as follows:

[0156]

[0157] NDVI was used as the greenness index, and LST was used as the heat index.

[0158] D. Create time series data for NDVI, NDWI, IBI, and SI respectively, and perform Mann-Kendall trend tests on each. Create an invariant layer. Taking NDVI as an example, create a layer. The calculation formula is as follows:

[0159] UC = Mann-Kendall{NDVI i}, i = the time series being monitored.

[0160]

[0161] Compare with the above formula and establish respectively And create an invariant layer. UC The formula is as follows:

[0162]

[0163] E. Based on Invariant Layers UC Time series data for indicators such as NDVI (Greenness Detection Index), Wet (Humidity), NDBSI (Natural Desiccation Index), and LST (Large Temperature Index) were masked. For each layer of the masked time phase, probability space statistics ranging from 2% to 98% were performed, with 2% as the minimum and 98% as the maximum. The original unmasked layers of NDVI, Wet, NDBSI, and LST were then normalized. Taking 2020 NDVI as an example, the normalized NDVI layer is as follows: The calculation formula is as follows:

[0164]

[0165] The dryness and heat indices are inversely normalized. Taking heat LST as an example, the normalization is as follows:

[0166]

[0167] F. By randomly distributing 1000 sample points within the study area (the number of sample points can be set according to the study area's scope and data resolution), and extracting four normalized indicators from multiple time phases year by year based on these sample points, the four normalized indicators from the time series are then combined into four new indicator columns. and Based on the above index list, PCA principal component analysis was performed to obtain the loading coefficients of the first principal component, labeled as a, b, c, and d, respectively, corresponding to... and Four indicator columns.

[0168] G. For each phase of the observed time series, calculate its initial RSEI0 value, as follows:

[0169]

[0170] Where y represents the monitoring time phase.

[0171] H. Acquired for all observation phases Statistical maximum value RSEI across multiple time phasesmax and minimum value RSEI min Then, normalize the result to obtain the final RSEI, as shown in the following formula:

[0172]

[0173] I. Apply the regions of the yearly cloud mask from step B to the RSEI generated in step H. y The cloud mask area is then assigned an empirical value of 0.8.

[0174] This invention has undergone extensive experimental verification, and the experimental results are quite satisfactory. Figure 2 This image shows the monitoring results of the remote sensing ecological index in the Guangdong-Hong Kong-Macao Greater Bay Area from 1990 to 2020, based on the technical solution of this invention. Experimental results show:

[0175] (1) The overall ecological condition of the Greater Bay Area is good, with almost no extremely poor areas;

[0176] (2) In 1990, the ecologically poor areas were located in Guangzhou and Shenzhen;

[0177] (3) In 2000, the number of ecologically poor areas in Guangzhou, Shenzhen and Dongguan all increased significantly;

[0178] (4) In 2010, the number of ecologically poor areas continued to increase in Guangzhou, Foshan and Dongguan, while the number in Shenzhen changed to a decreasing trend;

[0179] (5) In 2020, the poor areas in Shenzhen almost disappeared, and the poor areas in Guangzhou, Foshan and Dongguan also showed a decreasing trend.

[0180] The monitoring results are basically consistent with the changes in the ecological status of the Guangdong-Hong Kong-Macao Greater Bay Area.

[0181] Example 2

[0182] According to another embodiment of the present invention, a remote sensing ecological index monitoring device suitable for long-term time series is provided, see [link to relevant documentation]. Figure 3 ,include:

[0183] The data synthesis unit 100 is used to synthesize the acquired long-term Landsat remote sensing data.

[0184] Remote sensing index calculation unit 200 is used to calculate the time series remote sensing index of Landsat remote sensing data after data synthesis;

[0185] Trend test unit 300 is used to perform Mann-Kendall trend test on the time series remote sensing indices of Landsat remote sensing data to detect invariant regions.

[0186] Normalization unit 400 is used to normalize the time series remote sensing index using the probability space of the invariant region;

[0187] Principal component analysis unit 500 is used to perform principal component analysis on normalized indices from multiple time phases to obtain unified principal component loading coefficients.

[0188] The remote sensing ecological index monitoring device applicable to long-term time series in this embodiment of the invention is based on medium spatial resolution Landsat remote sensing imagery. It uses methods such as time series data synthesis, time series invariant region detection, and time series principal component analysis to design a remote sensing ecological index applicable to regions with drastic changes over a long period of time. Ultimately, it accurately reflects the distribution and evolution of the ecological environment over a large area and a long time series. After extensive experimental verification, the remote sensing ecological index monitoring results of this invention have high accuracy, the operation process is simple and feasible, and it has considerable practical value.

[0189] The following detailed description, using specific embodiments, illustrates the remote sensing ecological index monitoring device of the present invention applicable to long-term data series:

[0190] To address the shortcomings of existing technologies, this invention tackles issues in current remote sensing ecological index monitoring of ecological conditions, such as single-temporal normalization and multi-temporal comparison. It utilizes long-term Landsat remote sensing data as the data source to establish a remote sensing ecological index monitoring method suitable for long-term, drastically changing regions. This invention uses regional vegetation growth windows for data synthesis and employs the Mann-Kendall trend test for time-series remote sensing indices (NDVI, NDWI, IBI, SI) to detect invariant areas. NDVI, humidity, IBI, SI, and LST are normalized using a probability space of 2%–98% for invariant areas to avoid distortion caused by extreme maximum or minimum values. In the principal component analysis step, the normalized indices from multiple temporal phases are uniformly analyzed to obtain unified principal component loading coefficients. For water bodies, water body detection is used to mask the water body, and in the final RSEI results, gaps in the water body data are filled by direct assignment.

[0191] This invention is applicable to the overall evaluation of regional ecological conditions in cloudy and rainy areas, covering a wide range and over a long period of time. For areas where the ecological conditions have changed drastically, it can accurately express both the relative characteristics of the regional ecological conditions in a single time phase and the temporal changes of the ecological conditions in multiple time phases. The comprehensive mapping results of multiple time phases can help analyze the spatiotemporal evolution of the regional ecological conditions over a long period of time.

[0192] Specifically, this invention, based on medium spatial resolution Landsat remote sensing imagery, utilizes methods such as time-series data synthesis, time-series invariant region detection, and time-series principal component analysis to design a remote sensing ecological index suitable for regions experiencing drastic long-term changes. Ultimately, it accurately reflects the distribution and evolution of the ecological environment over a large area and over a long time period. Extensive experimental verification has demonstrated that the remote sensing ecological index of this invention achieves high accuracy, has a simple and feasible operation process, and possesses considerable practical value. The detailed technical solution of this invention is described below:

[0193] A. Data Synthesis Unit 100: Inputs time-series Landsat image data of the same region, performs cloud cover detection on each image, and masks the cloud cover pixels. Selects the peak vegetation growth period of the past 3 years for data synthesis. Taking South China as an example, the peak growth period is selected as September-October. The calculation formula is as follows:

[0194]

[0195] Where y represents the year of synthesis, and y-1≤y i ≤y+1, d1≤d j ≤d2, where d1 represents September 1st, d2 represents October 31st, and median is the median operator.

[0196] B. Perform water body detection on each synthesized Landsat image, mask the water-covered pixels, and save the water mask file for each image.

[0197] C. Remote Sensing Index Calculation Unit 200: For each composite Landsat image, calculate the remote sensing indices: NDVI, NDWI, IBI, SI, and temperature LST. The calculation formulas are as follows:

[0198] NDVI = (B nir -B red ) / (B nir +B red )

[0199] NDWI = (B green -B nir ) / (B green +B nir )

[0200]

[0201]

[0202]

[0203] in L6 = gain × DN + bias;

[0204] L6 is the radiometrically calibrated thermal infrared reflectivity. For the TM sensor, K1 = 607.76 W·m⁻²·sr⁻¹·μm⁻¹, K2 = 1260.56 K; for the TIRS sensor in band 10, K1 = 774.89 W·m⁻²·sr⁻¹·μm⁻¹, K2 = 1321.08 K; λ is the center wavelength of the thermal infrared band. TM =11.435μm, λ b10 =10.900μm; ρ =1.438×10-2m·K; ε is the surface emissivity.

[0205] Based on the tassel transformation, the humidity index Wet is established as follows:

[0206] Wet TM =0.0315×B blue +0.2021×B green +0.3012×Bρ red +0.1594×B nir -0.6806×B swir1 -0.6109B swir2

[0207] Wet OLI =0.1511×B blue +0.1973×B green +0.3283×B red +0.3407×B nir -0.7117×B swir1 -0.4559B swir2

[0208] Wet TM Wet OLI These correspond to the Landsat5TM sensor and the Landsat8OLI sensor, respectively.

[0209] Based on IBI and SI, the dryness index NDBSI is established, and the formula is as follows:

[0210]

[0211] NDVI was used as the greenness index, and LST was used as the heat index.

[0212] D. Trend Test Unit 300: Mann-Kendall trend tests are performed on NDVI, NDWI, IBI, and SI time series data respectively. An invariant layer is created; taking NDVI as an example, a layer is created. The calculation formula is as follows:

[0213] UC = Mann-Kendall{NDVI i}, i = the time series being monitored.

[0214]

[0215] Compare with the above formula and establish respectively And create an invariant layer. UC The formula is as follows:

[0216]

[0217] E. Normalized Unit 400: Based on Invariant Layer UC Time series data for indicators such as NDVI (Greenness Detection Index), Wet (Humidity), NDBSI (Natural Desiccation Index), and LST (Large Temperature Index) were masked. For each layer of the masked time phase, probability space statistics ranging from 2% to 98% were performed, with 2% as the minimum and 98% as the maximum. The original unmasked layers of NDVI, Wet, NDBSI, and LST were then normalized. Taking 2020 NDVI as an example, the normalized NDVI layer is as follows: The calculation formula is as follows:

[0218]

[0219] The dryness and heat indices are inversely normalized. Taking heat LST as an example, the normalization is as follows:

[0220]

[0221] F. Principal Component Analysis Unit 500: By randomly distributing 1000 sample points within the study area (the number of sample points can be set according to the study area range and data resolution), four normalized indicators from multiple time phases are extracted from the sample points year by year. These four normalized indicators from the time series are then combined into four new indicator columns. and Based on the above index list, PCA principal component analysis was performed to obtain the loading coefficients of the first principal component, labeled as a, b, c, and d, respectively, corresponding to... and Four indicator columns.

[0222] G. For each phase of the observed time series, calculate its initial RSEI0 value, as follows:

[0223]

[0224] Where y represents the monitoring time phase.

[0225] H. Acquired for all observation phases Statistical maximum value RSEI across multiple time phases max and minimum value RSEI min Then, normalize the result to obtain the final RSEI, as shown in the following formula:

[0226]

[0227] I. Apply the regions of the yearly cloud mask from step B to the RSEI generated in step H. y The cloud mask area is then assigned an empirical value of 0.8.

[0228] This invention has undergone extensive experimental verification, and the experimental results are quite satisfactory. Figure 2 This image shows the monitoring results of the remote sensing ecological index in the Guangdong-Hong Kong-Macao Greater Bay Area from 1990 to 2020, based on the technical solution of this invention. Experimental results show:

[0229] (1) The overall ecological condition of the Greater Bay Area is good, with almost no extremely poor areas;

[0230] (2) In 1990, the ecologically poor areas were located in Guangzhou and Shenzhen;

[0231] (3) In 2000, the number of ecologically poor areas in Guangzhou, Shenzhen and Dongguan all increased significantly;

[0232] (4) In 2010, the number of ecologically poor areas continued to increase in Guangzhou, Foshan and Dongguan, while the number in Shenzhen changed to a decreasing trend;

[0233] (5) In 2020, the poor areas in Shenzhen almost disappeared, and the poor areas in Guangzhou, Foshan and Dongguan also showed a decreasing trend.

[0234] The monitoring results are basically consistent with the changes in the ecological status of the Guangdong-Hong Kong-Macao Greater Bay Area.

[0235] Example 3

[0236] A storage medium storing program files capable of implementing any of the above-mentioned remote sensing ecological index monitoring methods applicable to long-term time series.

[0237] Example 4

[0238] A processor for running a program, wherein the program executes any of the above-mentioned remote sensing ecological index monitoring methods applicable to long-term series.

[0239] To eliminate discrepancies between RSEI values ​​and actual ecological conditions in water areas monitored in a single time phase, this invention extracts and masks water bodies from images of each time phase. The RSEI values ​​of the masked water areas are then filled in with empirical values ​​in the final step. To avoid normalization distortion caused by extreme points, this invention uses 2% and 98% values ​​from the statistical probability space to represent the minimum and maximum values, respectively, for index normalization. To eliminate errors in comparing multiple time phases in areas of drastic ecological change, invariant regions are detected and acquired from data of continuous time phases. A unified principal component analysis result is generated based on the invariant region indexes from multiple time phases and applied to all observation time phases to ensure consistency between the time series changes and the actual situation.

[0240] It is worth noting that:

[0241] (1) The original method did not process the water body. In this invention, the water body is extracted and masked in advance, and the extraction and masking are performed separately for each time phase. After the final RSEI is generated, the water body mask part is filled with empirical values ​​to ensure the integrity of the data spatial distribution.

[0242] (2) In the original method, the maximum and minimum values ​​are used to normalize each index. In this invention, the 2% and 98% values ​​of the statistical probability space are used to represent the minimum and maximum values, respectively. Moreover, reverse normalization is used for the dryness and heat indexes.

[0243] (3) In the original method, each time phase is normalized based on all its statistical characteristics. In this invention, on the one hand, the region participating in the feature statistics is determined by extracting the invariant region. On the other hand, a unified PCA principal component analysis is performed by random sampling and integration of all layers to generate a unified first principal component loading coefficient, which is then applied to the final RSEI calculation. This minimizes the impact of drastic changes in the long-term series region and can more accurately reflect the changes in the time series.

[0244] The sequence numbers of the above embodiments of the present invention are for descriptive purposes only and do not represent the superiority or inferiority of the embodiments.

[0245] In the above embodiments of the present invention, the descriptions of each embodiment have different focuses. For parts not described in detail in a certain embodiment, please refer to the relevant descriptions of other embodiments.

[0246] In the several embodiments provided in this application, it should be understood that the disclosed technical content can be implemented in other ways. The system embodiments described above are merely illustrative; for example, the division of units can be a logical functional division, and in actual implementation, there may be other division methods. For example, multiple units or components may be combined or integrated into another system, or some features may be ignored or not executed. Furthermore, the coupling or direct coupling or communication connection shown or discussed may be through some interfaces; the indirect coupling or communication connection of units or modules may be electrical or other forms.

[0247] The units described as separate components may or may not be physically separate. The components shown as units may or may not be physical units; that is, they may be located in one place or distributed across multiple units. Some or all of the units can be selected to achieve the purpose of this embodiment according to actual needs.

[0248] Furthermore, the functional units in the various embodiments of the present invention can be integrated into one processing unit, or each unit can exist physically separately, or two or more units can be integrated into one unit. The integrated unit can be implemented in hardware or as a software functional unit.

[0249] If the integrated unit is implemented as a software functional unit and sold or used as an independent product, it can be stored in a computer-readable storage medium. Based on this understanding, the technical solution of the present invention, in essence, or the part that contributes to the prior art, or all or part of the technical solution, can be embodied in the form of a software product. This computer software product is stored in a storage medium and includes several instructions to cause a computer device (which may be a personal computer, server, or network device, etc.) to execute all or part of the steps of the methods of the various embodiments of the present invention. The aforementioned storage medium includes various media capable of storing program code, such as USB flash drives, read-only memory (ROM), random access memory (RAM), portable hard drives, magnetic disks, or optical disks.

[0250] The above description is only a preferred embodiment of the present invention. It should be noted that for those skilled in the art, several improvements and modifications can be made without departing from the principle of the present invention, and these improvements and modifications should also be considered within the scope of protection of the present invention.

Claims

1. A remote sensing ecological index monitoring method suitable for long-term time series, characterized in that, Includes the following steps: Data synthesis was performed on the acquired long-term Landsat remote sensing data; Calculate the time-series remote sensing index of the Landsat remote sensing data after data synthesis; Mann-Kendall trend test was performed on the time series remote sensing indices of Landsat remote sensing data to detect invariant regions. Normalization of time series remote sensing indices is performed using the probability space of invariant regions; By uniformly performing principal component analysis on normalized indices from multiple time phases, a unified principal component loading coefficient is obtained; where: The method of performing the Mann-Kendall trend test on the time-series remote sensing indices of Landsat remote sensing data to detect invariant regions includes: Mann-Kendall trend tests were performed on NDVI, NDWI, IBI, and SI time series data respectively, and invariant layers were created. A layer was created based on the NDVI time series data. The calculation formula is as follows: Compare with the above formula and establish respectively , , And create an invariant layer. The formula is as follows: ； The normalization of the time-series remote sensing index using the probability space of the invariant region includes: Based on invariant layers Time series data for NDVI (Greenness Spectrum Indices), Wet (Humidity Spectrum Indices), NDBSI (Dryness Spectrum Indices), and LST (Heat Spectrum Indices) were masked. For each layer of the masked time phase, probability space statistics ranging from 2% to 98% were performed, with 2% as the minimum and 98% as the maximum. The original unmasked layers for greenness, humidity, dryness, and heat were then normalized. The normalized greenness layer for NDVI in 2020 was [details omitted]. The calculation formula is as follows: The dryness and calorific value indices are inversely normalized, with the calorific value LST being inversely normalized. The calculation formula is as follows: ； The step of performing principal component analysis on normalized indices from multiple time phases to obtain unified principal component loading coefficients includes: Several sample points were randomly deployed within the study area, with the number of sample points determined based on the study area's extent and data resolution. Based on these sample points, four normalized indicators from multiple time phases were extracted annually. These four normalized indicators from the time series were then combined into four new indicator columns. , , and Based on the above index list, PCA principal component analysis was performed to obtain the loading coefficients of the first principal component, labeled as a, b, c, and d, respectively, corresponding to... , , and 4 indicator columns; After performing principal component analysis on the normalized indices of multiple time phases to obtain unified principal component loading coefficients, the method further includes: For each phase of the observed time series, the initial RSEI0 value is calculated using the following formula: Where y represents the monitoring time phase; Acquired for all observation phases Statistical maximum value of multiple time phases and minimum value Then, normalize the result to obtain the final RSEI, as shown in the following formula: 。 2. The remote sensing ecological index monitoring method applicable to long-term time series according to claim 1, characterized in that, The data synthesis of the acquired long-term Landsat remote sensing data includes: Input Landsat time-series image data of the same region, perform cloud cover detection on each image, mask the cloud cover pixels, and synthesize the data by selecting the nearest N years of peak vegetation growth period. The calculation formula is as follows: in, Indicates the year of synthesis. , , It is the median operator.

3. The remote sensing ecological index monitoring method applicable to long-term time series according to claim 2, characterized in that, After synthesizing the acquired long-term Landsat remote sensing data, the method further includes: For each synthesized Landsat image, water bodies are detected, and the pixels covered by water bodies are masked.

4. The remote sensing ecological index monitoring method applicable to long-term time series according to claim 3, characterized in that, The time-series remote sensing indices of the synthesized Landsat remote sensing data include: For each composite Landsat image, the remote sensing indices NDVI, NDWI, IBI, SI, and temperature LST were calculated using the following formulas: in , ; K1 represents the reflectivity of the thermal infrared band after radiometric calibration. For the TM sensor, K1 = 607.76 W·m⁻²·sr⁻¹·μm⁻¹, K2 = 1260.56 K; for the TIRS sensor in band 10, K1 = 774.89 W·m⁻²·sr⁻¹·μm⁻¹, K2 = 1321.08 K; λ is the center wavelength of the thermal infrared band. , ; 10⁻² m·K; The surface emissivity; Based on the tassel transformation, a humidity index Wet is established, with the following formula: , These correspond to the Landsat5TM sensor and the Landsat8OLI sensor, respectively. Based on IBI and SI, the dryness index NDBSI is established, and the formula is as follows: NDVI was used as the greenness index, and LST was used as the heat index.

5. The remote sensing ecological index monitoring method applicable to long-term time series according to claim 4, characterized in that, After performing principal component analysis on the normalized indices of multiple time phases to obtain unified principal component loading coefficients, the method further includes: The regions of the yearly cloud masks, after masking the water-covered pixels, were applied to all observation phases. Statistical maximum value of multiple time phases and minimum value Then, normalize the result to obtain the final RSEI generated in the process. And assign the cloud mask area an empirical value.

6. A remote sensing ecological index monitoring device suitable for long-term series, utilizing the remote sensing ecological index monitoring method for long-term series as described in claim 1, characterized in that, include: The data synthesis unit is used to synthesize the acquired long-term Landsat remote sensing data. The remote sensing index calculation unit is used to calculate the time series remote sensing index of the Landsat remote sensing data after data synthesis. The trend test unit is used to perform the Mann-Kendall trend test on the time series remote sensing indices of Landsat remote sensing data to detect invariant regions. The normalization unit is used to normalize the time series remote sensing index using the probability space of the invariant region; The principal component analysis unit is used to perform principal component analysis on normalized indices from multiple time phases to obtain unified principal component loading coefficients.

Images