Weakly supervised temporal feature driven remote sensing cropland product fusion mapping method

Abstract

The application discloses a weakly supervised time sequence feature driven remote sensing cultivated land product fusion mapping method, comprising the following steps: weakly supervised time sequence dataset construction; constructing a utae time sequence network model which is reformed by a weakly supervised loss; weakly supervised model training, after obtaining a preliminary model, supplementary training is carried out on difficult samples; finally, the trained model is used for cultivated land coverage mapping in a proper range. The weakly supervised learning idea is applied to the cultivated land product fusion mapping task, and the product fusion precision is improved by introducing the time sequence model.

Description

TECHNICAL FIELD

[0001] The application belongs to the technical field of remote sensing cultivated land product mapping, and particularly relates to a weakly supervised time sequence feature driven remote sensing cultivated land product fusion mapping method. BACKGROUND

[0002] As a basic resource for human survival and social development, the area and distribution of cultivated land are closely monitored by the Chinese government, and the red line of 1.8 billion mu of cultivated land must be strictly guarded. Timely and accurate cultivated land mapping is crucial. With the rapid development of remote sensing technology, a large number of satellite sensors with different spatial and temporal resolutions and spectral resolutions have been launched, providing opportunities for faster and more accurate monitoring of surface resources. Currently, there are two types of products that can serve cultivated land monitoring: cultivated land thematic mapping products and global land cover maps.

[0003] For cultivated land thematic mapping, countries have launched a series of related projects. Common cultivated land mapping products include the CDL (Cropland Data Layer) crop classification map in the United States, the AFCC product in Canada, and the Sen-Agri in Europe. Since the production of these maps often requires sufficient labeled data for model training and a large number of field surveys to verify product reliability, the coverage of the products is limited to the implementation areas of the above projects, and it is usually difficult to produce cultivated land thematic maps on a large scale.

[0004] For global land cover products, many regional or global land cover products have been released by different countries and organizations. For example, Google's 10m near real-time global land cover project: Dynamic World, Esri's 2017-2021 five-year land cover product, Tsinghua University's FROM-GLC, etc. The emergence of these high spatial resolution land cover products provides an opportunity for global-scale cultivated land monitoring. However, to ensure global coverage of the product and overall accuracy of all land cover types, multiple images are often used to splice single temporal data or data with long time intervals for classification. The optical characteristics of cultivated land change significantly with the seasons, and the visual features such as texture and color of cultivated land images in the same area at different times also change greatly. Therefore, the accuracy of cultivated land in the land cover results obtained by mapping using single temporal remote sensing images is not high.

[0005] Because the distribution and area of cultivated land generally do not change dramatically over time, in order to produce and update global cultivated land coverage products with high precision at low cost, many researchers have attempted to fuse multiple land cover products according to certain standards to generate new cultivated land thematic maps with higher precision. This paper introduces the current status and technical deficiencies of land cover product fusion technology from two aspects: data sources and product fusion methods.

[0006] Product fusion data sources: The spatial resolution of remote sensing images limits the mapping precision of land cover products. Early land cover products were mainly at low spatial resolution. With the launch of the European Sentinel series of satellites and the opening of shared data from China's Gaofen and Resource series of satellites, land cover product production has entered the medium and high spatial resolution era. With the development of machine learning methods and remote sensing cloud computing platforms, high spatial resolution and high frequency global land cover mapping has become possible. The current commonly used global land cover products are as follows:

[0007] Google's Dynamic World product: Google released this product in 2022, with a spatial resolution of 10 m, covering the global range, and a time span of 2015-2022. The product uses Sentinel 2 Level 1C products as images. It uses a deep learning model as the classification algorithm, with 24,000 globally labeled images as the training basis, and relies on Google's own cloud platform with high computing power to achieve near real-time land cover mapping.

[0008] ESA's World Cover product: ESA released this product in 2021. The WC product releases two periods in 2020 and 2021, with a spatial resolution of 10 m, covering the global range. WC uses Sentinel 2 and Sentinel 1 as the basis for data, and uses open source geographic data such as OSM as auxiliary data. It uses a random forest classification tree algorithm to train and predict on manually labeled pixels in 100x100 meter grids.

[0009] Esri's Land Cover product: Esri released this product in 2021, with a spatial resolution of 10 m, a time span of 2017-2021, and a coverage of the global range. The product uses Sentinel 2 as the model input data, and each period uses single temporal data. Esri's product uses deep learning methods for production, and uses the same training data set as Google's product, but after classification, it fine-tunes the classification results according to priori.

[0010] However, on the one hand, the production process of these products does not take into account the spatio-temporal heterogeneity of cultivated land, on the other hand, due to the different land cover classification systems used between different data, different classification methods, different verification methods, there will be great inconsistency between different products when they are fused, and direct superposition of different products may bring different product error classification results to the fusion product, thus reducing the fusion accuracy.

[0011] Product fusion technology method: product fusion method needs to overcome the defects of existing land cover products and the incompatibility between different products, in order to get higher precision fusion results. Current product fusion methods can be classified into learning and non-learning product fusion methods.

[0012] Non-learning: Non-learning method tries to combine statistical data fusion method to solve the problem of compatibility, comparability and accuracy of existing products. According to the fusion method, these studies can be divided into GWR-based and fusion decision rule-based. The feature of GWR method is that the regression parameter changes with the change of spatial distance, so it can better reflect the influence of spatial position factors. Linda See et al. in the non-patent literature "Building a hybrid land cover map with crowdsourcing and geographically weighted regression, ISPRS Journal of Photogrammetry and Remote Sensing, 103:48-56, 2015" fused GLC2000, MODIS global land cover product, GlobCover three products by GWR method, and generated a global 300m product map.

[0013] Fusion decision rules include various models, Bayesian theory, Dempster-Shafer evidence theory. Bayesian theory can improve the separability between classes under the condition of known prior probability distribution. A. Perez-Hoyos et al. used Bayesian theory to make a fused global land cover map in the non-patent document: Integrating multiple land cover maps through a multi-criteria analysis to improve agricultural monitoring in Africa, International Journal of Applied Earth Observation and Geoinformation, 2020, and gave the confidence value of each pixel classification, reducing the error caused by data coordination and spatial sampling and other factors.

[0014] Learning type: The learning type product fusion method attempts to collect different product samples with high reliability as labels into a classifier, and regenerate a more accurate local or global map through additional remote sensing images. Due to the uneven quality of products, a sample screening mechanism is usually needed to select reliable samples. Common methods include rule-based screening, such as only keeping samples that have not changed over multiple years, only keeping central pixels to reduce the influence of edge effects, or using edge reduction operators to remove boundary pixels. Liu et al. filtered pixels and spectra from MODIS land cover products through GEE (Google Earth Engine) in the non-patent document “High-resolution multi-temporal mapping of global urban land using Landsat images based on the Google Earth Engine Platform, Remote sensing of Environment, 2018, 209:227-239” to screen labels, and finally achieved a sample accuracy of 99.2%.

[0015] However, on the one hand, the non-learning data fusion method can reduce the problem of incoordination and incompatibility between different products, but lacks remote sensing image assistance and cannot be optimized for cultivated land thematic map. The result obtained by using the original product is still a single time sequence result. On the other hand, in order to reduce the error of the product itself and the product, the existing learning fusion framework needs to discard the inconsistent pixels between different products, so that complete and accurate labels required by deep learning cannot be obtained. The supervision classification method of deep learning usually needs complete and accurate labels as model guidance to train a reliable classification model. This also leads to the fact that the traditional deep learning method cannot be directly used for learning product fusion method, and the current learning method is mainly based on machine learning. Finally, neither the non-learning data fusion framework nor the learning data fusion framework considers the time sequence characteristics of cultivated land. Single time sequence image data is still used in the fusion process, and the final accuracy cannot be further improved. SUMMARY

[0016] Therefore, the purpose of the present application is to fuse the current land cover product with higher accuracy to quickly produce cultivated land cover product. To this end, on the one hand, the weak supervision learning method is introduced into the current product fusion framework to improve the learning product fusion method, so as to learn a reliable model for cultivated land cover mapping in the case of incomplete labels. On the other hand, long time sequence image data is used as the remote sensing data source for cultivated land time sequence feature modeling. A utae network, a deep learning network with time sequence feature coding structure, is used as the basic model to model the time sequence features of cultivated land to improve the final cultivated land cover mapping accuracy. Based on the above two points, the present application proposes a remote sensing cultivated land product fusion mapping framework based on weak supervision learning time sequence feature extraction.

[0017] Specifically, the weak supervision time sequence feature driven remote sensing cultivated land product fusion mapping method proposed by the present application comprises the following steps:

[0018] Weak supervision time sequence dataset construction;

[0019] Constructing a utae time sequence network model modified by a weak supervision loss;

[0020] Weak supervision model training, after obtaining the preliminary model, the difficult samples are further trained;

[0021] Finally, the trained model is used for cultivated land cover mapping in a proper range.

[0022] Further, the utae time sequence network model modified by the weak supervision loss comprises:

[0023] Step 201: Obtain the output of the last layer of the utae decoding layer, that is, the decoded image feature vector V before the network softmax layer f At this time, the size of the feature vector is C*W*H, where W and H are spatial dimensions, and C is a channel dimension; In the spatial dimension, a feature pixel set S = {p1, p2, p3...p n} with a size of 0.01*W*H is randomly selected at a proportion of 1%; The size of each feature is C*1*1;

[0024] Step 202: Calculate the spatial dimension correlation: for the feature set S obtained in step 201, each element is operated as follows:

[0025] Take the feature p i The 8 pixels adjacent in the spatial dimension, if on the image boundary, only 3-5 pixels; Calculate the Euclidean distance between the element and the feature of the adjacent pixel in the spatial dimension, and select the two results closest and farthest as the two feature vectors closest and farthest to the pixel;

[0026] Step 203: Calculate the feature dimension correlation: for the feature set S obtained in step 201, each element is operated as follows:

[0027] Take the feature p i A new channel feature set S v is formed by multiple feature vectors V' in the channel dimension, where the size of each vector is 1*1*1, and the cosine similarity between each element in S v is calculated in turn, and the minimum value is taken as the minimum similarity C min in the channel feature;

[0028] Step 204: Construct weak supervision loss: for each element in the feature set S, the minimum and maximum spatial similarity distances d min and d max obtained by the calculation in step 202 and step 203, and the minimum feature similarity C min obtained in the channel dimension are summed to obtain the weak supervision loss:

[0029]

[0030] Where α, β, γ are weight parameters of each feature similarity, used to adjust the contribution of each similarity to the final loss;

[0031] Step 205: Construct the complete loss of the final model: only calculate the part of the label value 0 and 1, ignore the pixels with the label grid value 2, and ignore the pixels with the label grid value 2, and the cross entropy calculation formula is:

[0032] loss ce =Σ-(y i ·log(p i )+(1-y i )·log(1-p i ))

[0033]

[0034] Where y i is the label corresponding grid value, is the network prediction result probability;

[0035] The final loss function formula is as follows:

[0036] loss=loss ce +λ·loss weak

[0037] Where λ is used to weight the weakly supervised loss, and the weakly supervised loss calculation is based on the feature vector of the model after the decoding layer.

[0038] Further, the discrete weakly supervised label and the sentinel-2 time series data are used to supervise the learning of the reconstructed utae time series model;

[0039] The basic utae network adds an L-TAE module based on the u-net network to perform time series feature encoding. The L-TAE is a simplified multi-head self-attention mechanism network that directly processes the deep abstract features of the image after convolution encoding, thereby transforming the basic u-net model into a time series utae model that can process time series.

[0040] Further, the time series feature encoding step is as follows:

[0041] Step 301: Perform L-layer convolution operation on the time series data to obtain the deep feature map e L of the original image, and the size of e L is T*C L *W L *H L :

[0042] e l =ε l (e l-1 )for l∈[1,L]

[0043] Where ε l() represents the convolution operation on the characteristics of the l-1 layer, and the depth of the convolution characteristics is L when l=L;

[0044] Step 302: The L-TAE multi-head self-attention module acts on the feature map e L Temporal feature encoding is performed to obtain G groups of encoded temporal feature weights Each group of weights has a size of T*W L *H L :

[0045]

[0046] Step 303: In order to use the temporal feature weights on different depth encoders, the temporal feature weights calculated in step 302 are resampled to have the same length and width as the shallow layer features:

[0047] a l,g =resize l (a L,g )

[0048] Where resize l () represents resampling the feature in the spatial dimension to W l *H l , and a l,g is the resampled temporal feature weight;

[0049] Step 304: The temporal feature weight a l,g is applied to the convolution feature map e l obtained by each layer of the encoder, and the action process is as follows:

[0050] Each layer of features is evenly grouped to obtain G groups of feature vectors with a size of T*C l / G*W l *H l The feature vectors are multiplied by the temporal weight a l,g term by term, and summed in the time dimension, and a convolution operation is performed on the result to obtain the temporally encoded feature map:

[0051]

[0052] Further, the model is used for cultivated land coverage mapping. In the mapping link, the cultivated land coverage is predicted using the temporally incomplete image data. The 256x256 size samples obtained by sliding cropping are used as units for prediction during mapping. The part with higher consistency of cultivated land class in the three products, i.e. the part with a grid value of 1 in the weakly supervised data set, is retained. The result predicted by the weakly supervised model is used to modify the inconsistent area of the product, and finally the complete cultivated land coverage mapping result is obtained by splicing.​

[0053] The beneficial effects of the present application are as follows:

[0054] (1) The present application integrates the idea of weak supervision into the framework of cultivated land product fusion mapping, solving the problem of a large number of effective pixel losses caused by inconsistency between products in the existing product fusion framework. And in view of the existence of a large number of existing free products and free remote sensing images, therefore, the method can be theoretically promoted in the global cultivated land mapping.

[0055] (2) For the problem of cultivated land thematic mapping, considering the seasonal change of cultivated land over time, the present application proposes a cultivated land product fusion mapping framework based on time series image data, which uses sentinel 2 time series data as the image data source, and uses the weak supervision module to modify the utae model to encode the time series characteristics of the cultivated land image, so that the mapping result is more in line with the time series characteristics of the cultivated land, and effectively improves the accuracy of the cultivated land type in the previous land cover product. BRIEF DESCRIPTION OF DRAWINGS

[0056] Figure 1 Figure 1 is the flowchart of the present application;

[0057] Figure 2 is the product fusion framework diagram proposed by the present application;

[0058] Figure 3 is a weak supervision label superposition diagram;

[0059] Figure 4 is a weak supervision learning module modified utae model used by the present application. DETAILED DESCRIPTION

[0060] The present application will be further described below in conjunction with the drawings, but in any way limits the present application, any transformation or replacement based on the teaching of the present application belongs to the protection scope of the present application.

[0061] To achieve the above-mentioned purposes, the present application adopts the following technical solutions:

[0062] Step 1, weak supervision time series dataset production;

[0063] Step 2, construct the utae time series network model modified by the weak supervision loss;

[0064] Step 3, weak supervision model training and difficult sample supplementary training;

[0065] Step 4, using the trained model to carry out the cultivated land cover mapping in the appropriate range.

[0066] Weakly supervised learning can include the following three modes: incomplete label training: only part of the samples in the training set are labeled, and most of the samples are not labeled; imprecise label training: the samples in the training set only have rough labels, such as semantic segmentation labels using only line labels for annotation; inaccurate label training: part of the samples in the training set are labeled incorrectly. The weakly supervised label referred to in the present application is an imprecise label, that is, only part of the cultivated land in a sample is labeled.

[0067] Weakly supervised time series dataset: according to the needs of the study of cultivated land map, the corresponding region of the time series image and multiple land cover products in the geographical spatial range are downloaded, and in the case of the present application, Sentinel-2 multispectral data is used as the time series image, Google's DynamicWorld product, ESA's WorldCover product, and Esri's global land cover product are used as the basic land cover data source. Finally, a weakly supervised dataset of incomplete cultivated land semantic segmentation labels with an image time series length of 12 months is obtained.

[0068] Constructing a utae time series network modified by a weakly supervised loss, the original utae network needs to provide complete and accurate semantic segmentation labels, so it cannot be directly used with the weakly supervised labels described in the present application. The output part of the network and the loss function need to be modified to learn more accurate results on the subsequent weakly supervised labels. The specific steps are as follows:

[0069] Step 201: Obtain the output of the last layer of the utae decoding layer, that is, the decoded image feature vector V before the softmax layer of the network f At this time, the size of the feature vector is C*W*H, where W and H are spatial dimensions, and C is the channel dimension. In the spatial dimension, 1% of the feature pixels are randomly selected to obtain a feature set S = {p1, p2, p3...p n}n = 0.01 x W x H, where the size of each feature is C*1*1.

[0070] Step 202: Calculate the spatial dimension correlation: for the feature set S obtained in step 201, perform the following operations on each element in it:

[0071] Take the feature P i The 8 pixels P i L ,P i R ,P i T ,P i B ,P i LT ,Pi RT P i RB P i LB If on the image border, only 3-5 pixels. Calculate the Euclidean distance of this element with its spatial neighbors and select the two results with the closest and farthest distances:

[0072] d max = Max(D(P i , P i L ), D(P i , P i R ), D(P i , P i T ), D(P i , P i B ), D(P i , P i LT ), D(P i , P i RT ), D(P i , P i RB ), D(P i , P i LB ))

[0073] d min = Min(D(P i , P i L ), D(P i , P i R ), D(P i , P i T ), D(P i , P i B ), D(P i , P i LT ), D(P i , P i RT ), D(P i , P i RB ), D(P i , P i LB ))

[0074] As the two feature vectors P that are closest and farthest from this pixel min ,P max The formula for calculating Euclidean distance is as follows:

[0075] D = ||P i ,P j || 2

[0076] Step 203: Calculate feature dimension correlation: For the feature set S obtained in step 201, perform the following operation on each element:

[0077] Take this feature P i Multiple feature vectors V' in the channel dimension constitute a new channel feature set S. v Each vector has a size of 1*1*1, and S is calculated sequentially. v The minimum cosine similarity between each element in the algorithm is taken as the minimum similarity C on the channel features. min The formula for calculating cosine similarity is:

[0078]

[0079] Where V i 'and V j 'represents from set S v The element extracted from the array.

[0080] Step 204: Construct the weakly supervised loss: For each element in the feature set S, calculate the minimum and maximum spatial similarity distance d between each element obtained in steps 202 and 203. min and d max The minimum feature similarity C calculated on the channel dimension min Summing the three eigenvalues ​​yields the weakly supervised loss:

[0081]

[0082] Where α, β, and γ are the weight parameters for each feature similarity, which are used to adjust the contribution of each similarity to the final loss.

[0083] Step 205: Construct the complete loss for the final model: Modify the general cross-entropy loss, stipulating that it is only calculated for the portion of the label values ​​0 and 1, ignoring pixels with a raster value of 2 on the label. Cross-entropy calculation formula:

[0084] loss ce =Σ-(y i ·log(p i )+(1-y i )·log(1-pi ))

[0085] Where y i p is the raster value corresponding to the label. i This represents the probability of the network prediction result.

[0086] The final loss function formula is as follows:

[0087] loss = loss ce +λ·loss weak

[0088] Here, λ is used to weight the weakly supervised loss, and in the results obtained in this invention, λ is set to 0.1. It is worth noting that the calculation of cross-entropy loss is based on the final result of the model after passing through the softmax layer, while the calculation of weakly supervised loss is based on the feature vector of the model after passing through the decoding layer.

[0089] Discrete weakly supervised labels and Sentinel-2 time-series data are used to supervise the learning of the modified UTAE time-series model. The basic UTAE network is a typical U-shaped network structure. An L-TAE (Lightweight-Temporal Attention Encoder) module is added to the U-Net network for temporal feature encoding. L-TAE is a simplified multi-head self-attention mechanism network that can directly process the deep abstract features of images after convolutional encoding, thus transforming the basic U-Net model into a time-series UTAE model capable of processing time series data. (Reference: Lightweight temporal self-attention for classifying satellite images time series). The time encoding steps are as follows:

[0090] Step 301: Perform L layers of convolution on the time-series data to obtain the depth feature map e of the original image. L e L The size is T*C L *W L *H L :

[0091] e l =ε l (e l-1 for l∈[1,L]

[0092] Where ε l () indicates that a convolution operation is performed on the features of layer l-1. When l = L, a convolution feature of depth L is obtained.

[0093] Step 302: The L-TAE multi-head self-attention module processes the feature map e obtained in step 301.L Perform temporal feature encoding to obtain G sets of encoded temporal feature weights. The size of each weight group is T*W. L *H L :

[0094]

[0095] Step 303: In order to use temporal feature weights on encoders at different depths, the temporal feature weights calculated in step 302 are resampled to have the same width and height as the shallow layer features:

[0096] a l,g =resize l (a L,g )

[0097] where resize l () indicates that the features are resampled to W in the spatial dimension. l *H l .

[0098] Step 304: Weight the time series features a l,g The convolutional feature map e obtained by each encoder layer l The process is as follows:

[0099] Each layer of features is uniformly grouped to obtain G groups of size T*C. l / G*W l *H l eigenvectors Using time weight a l,g For eigenvectors Perform term-by-term multiplication and summation over time, then perform a convolution operation on the result to obtain the temporally encoded feature map:

[0100]

[0101] The parameters set during training are as follows:

[0102] The training epochs are 50, the batch size is set to 8, the initial learning rate is 0.001, the Adam optimizer is used, and the weight decay is 1e. -5 .

[0103] Additional training was performed on the model for difficult samples. Because some farmland samples had low recognition accuracy for two or more farmland products, the number of pixels for these types of farmland was extremely small during the weak supervision labeling process due to inconsistency screening. As a result, the final model did not recognize this type of farmland well. Therefore, about 30 samples were selected from the product that predicted this type of farmland well to allow the model to continue learning.

[0104] The model is used for farmland coverage mapping in a proper range. In the mapping link, the farmland coverage prediction can be performed using the time-series incomplete image data, and the 256x256 size sample obtained by sliding clipping is used as a unit for prediction. In the mapping process, the part with high consistency of the farmland class in the three products, i.e., the part with the grid value of 1 in the weak supervision data set, is retained, the result of the weak supervision model prediction is used to modify the inconsistent area of the product, and finally the complete farmland coverage mapping result is obtained by splicing.

[0105] The beneficial effects of the present application are as follows:

[0106] (1) The weak supervision idea is integrated into the farmland product fusion mapping framework, which solves the problem of a large number of effective pixel losses caused by inconsistency between products in the existing product fusion framework. And in view of the existence of a large number of existing free products and free remote sensing images, therefore, the method can be theoretically popularized in global farmland mapping.

[0107] (2) For the problem of farmland thematic mapping, considering the seasonal change of farmland over time, the present application proposes a farmland product fusion mapping framework based on time-series image data, which uses the sentinel 2 time-series data as the image data source, and uses the weak supervision module to modify the utae model to encode the time-series characteristics of the farmland image, so that the mapping result is more consistent with the time-series characteristics of the farmland, and effectively improves the accuracy of the farmland type in the previous land cover product.

[0108] Embodiment

[0109] As shown in Figure 1 , the weak supervision time-series feature driven remote sensing farmland product fusion mapping method comprises the following steps:

[0110] Step 1, weak supervision time-series data set making;

[0111] Step 2, constructing the utae time-series network model modified by the weak supervision loss;

[0112] Step 3, weak supervision model training and difficult sample supplementary training;

[0113] Step 4, using the trained model to perform farmland coverage mapping in a proper range.

[0114] In the present application, weakly supervised learning refers to using incomplete and sparse labeled labels in the semantic segmentation task to guide the model, and under the condition of incomplete labels, the consistency between the feature vectors of unknown pixels and known pixels is measured to determine whether the unknown pixels are related to the known pixels, so as to reduce the loss of supervised information caused by incomplete labels. In the process of land cover product fusion, due to different production methods of different products, different classification systems used in labeling, and different training sets used in model training, it is inevitable that there will be inconsistent classification results between different land cover products in some pixels. Such inconsistent pixels cannot be used as reliable pixels to learn in the model, but in fact, there are cases where one product is accurately predicted while the other product is incorrectly predicted. Therefore, such inconsistent pixels cannot be completely discarded as waste pixels. Therefore, the idea of weakly supervised learning is introduced into the product fusion framework, and the labels of the consistent area are used to guide the model learning in the learning process, and the inconsistent area is used as a weak supervision signal to indirectly affect the model. Considering that the cultivated land plot changes with the season and human activities, it has obvious different visual features at different times. In the fallow period, the cultivated land plot is close to bare soil, and in the mature period, it is close to the color of vegetation. In the harvest period, the green vegetation features of the cultivated land are quickly weakened under the influence of human activities, and in a short period of time, it becomes bare land. Therefore, it is difficult to accurately monitor the cultivated land range based on single temporal image data. In order to better fit the rapid change of cultivated land with the season, it is necessary to use time-series remote sensing image data and cooperate with the UTAE time-series model as the basic network to improve the prediction accuracy of cultivated land. Considering the inconsistency of the labels and the time-series characteristics of the cultivated land, a cultivated land product fusion framework as shown in Figure 2 is proposed. The framework details the entire process of the cultivated land product fusion described in the present application. The following describes the entire product fusion process in combination with the framework:

[0115] Weakly supervised time-series dataset making: According to the needs of making cultivated land map research, download the time-series images and multiple land cover products in the corresponding geographical space range. In the case of the present application, Sentinel-2 multispectral data is used as the time-series image, and Google's DynamicWorld product, ESA's WorldCover product, and Esri's global land cover product are used as the basic land cover data source. In the actual data collection process, even if the Sentinel-2 data has complete time-series data, there are still problems such as cloud cover in a single image, resulting in no data in a small area. This part of the data is displayed as null value in the grid data, so the mean value of the upper and lower two time-series images with null data is used to fill the small part of the null data area.

[0116] Regarding the making of weakly supervised labels, in combination with Figure 3The detailed process is as follows:

[0117] Step 101: Perform label mapping on the given three data source products P x ,P y ,P z , convert the land types in the product except the cultivated land into the same background pixel, retain the cultivated land pixel, and obtain the transformed cultivated land thematic product P' x ,P' y ,P' z . Perform raster addition on the three products to obtain a single raster product A with a raster value range of 0-3, perform label mapping on the raster A, map the pixel with a raster value of 1 to 2, and map the pixel with a raster value of 3 to 1 to obtain the final weakly supervised label L w , the part with a raster value of 2 in the above weakly supervised label does not participate in the training of the model, the part with a raster value of 1 is the cultivated land pixel, and 0 is the background pixel. The formulaic description of this step is as follows:

[0118] L w =f″(Sum(f′(P x ),f′(P y ),f′(P z )))

[0119] Where P x ,P y ,P z represent the original product, f' is the first label mapping rule, Sum represents raster addition, and f" is the second label mapping rule.

[0120] The utae time sequence network modified by the weakly supervised loss is constructed, the original utae network needs to provide complete and accurate semantic segmentation labels, therefore, the weakly supervised label described in the present application cannot be directly used, and the output part of the network and the loss function need to be modified to learn more accurate results on the subsequent weakly supervised label.

[0121] Referring to Figure 4 , the modification of the weakly supervised module is mainly in the output part of the network, and the specific steps are as follows:

[0122] Step 201: Obtain the output of the last layer of the utae decoding layer, that is, the decoded image feature vector V f before the softmax layer of the network, at this time, the feature vector size is C*W*H, wherein W and H are spatial dimensions, and C is a channel dimension. In the spatial dimension, feature pixels are randomly selected at a rate of 1%, and a feature vector set S = {p1, p2, p3...p n}n = 0.01*W*H is obtained, wherein the size of each feature is C*1*1.

[0123] Step 202: Compute spatial dimension correlation: For the feature set S obtained in step 201, for each element do the following:

[0124] Take the feature P i 8 pixels P i L ,P i R ,P i T ,P i B ,P i LT ,P i RT ,P i RB ,P i LB , if the pixel is on the image boundary, only 3-5 pixels. Compute the Euclidean distance between this element and its spatially adjacent pixel features, and select the two results with the closest and farthest distances:

[0125] d max = Max(D(P i ,P i L ), D(P i ,P i R ), D(P i ,P i T ), D(P i ,P i B ), D(P i ,P i LT ), D(P i ,P i RT ), D(P i ,P i RB ), D(P i ,P i LB ))

[0126] d min = Min(D(P i ,P i L ), D(P i ,P i R ), D(P i ,P i T ), D(Pi ,P i B ),D(P i ,P i LT ),D(P i ,P i RT ),D(P i ,P i RB ),D(P i ,P i LB ))

[0127] The two feature vectors P min ,P max closest and farthest to the pixel are taken as the minimum and maximum distance respectively.

[0128] The formula of the Euclidean distance is as follows:

[0129] D = ||P i ,P j || 2

[0130] Step 203: Calculate feature dimension correlation: for the feature set S obtained in step 201, the following operations are performed on each element in it:

[0131] Take the feature P i In the channel dimension, multiple feature vectors V' form a new channel feature set S v , where the size of each vector is 1*1*1, and the cosine similarity between each element in S v is calculated in turn, and the minimum value is taken as the minimum similarity C min on the channel feature. The formula for calculating the cosine similarity is:

[0132]

[0133] Where V i ' and V' j represent the elements taken from the set Sv.

[0134] Step 204: Construct weak supervision loss: for each element in the feature set S, the minimum and maximum spatial similarity distances d min and d max obtained in step 202 and step 203, and the minimum feature similarity C min obtained in the channel dimension are summed to obtain the weak supervision loss:

[0135]

[0136] where a, b, g are the weight parameters of each feature similarity, used to adjust the contribution of each similarity to the final loss.

[0137] Step 205: Construct the complete loss of the final model: modify the general cross-entropy loss, stipulate to calculate only in the part of label value 0 and 1, ignore the pixel whose label grid value is 2. The cross-entropy calculation formula is:

[0138] loss ce =∑-(y i ·log(p i )+(1-y i )·log(1-p i ))

[0139] where y i is the label corresponding to the grid value, and p i is the network prediction probability.

[0140] The final loss function formula is as follows:

[0141] loss=loss ce +λ·loss weak

[0142] Where λ is used to weight the weakly supervised loss, and the present application sets λ to 0.1. It is worth noting that the calculation of cross-entropy loss is based on the final result of the model after the softmax layer, while the calculation of weakly supervised loss is based on the feature vector of the model after the decoding layer.

[0143] The discrete weakly supervised label and sentinel 2 time series data are used to supervise the learning of the modified utae time series model. The basic utae network is a typical U-shaped structure network, which adds an L-TAE (Lightweight-Temporal Attention Encoder) module for time series feature encoding based on the u-net network. L-TAE is a simplified multi-head self-attention mechanism network that can directly process deep abstract features of images after convolution encoding, thereby modifying the basic u-net model into a time series utae model that can process time series. Reference: Lightweight temporal self-attention for classifying satellite images time series.

[0144] The time encoding step is as follows:

[0145] Step 301: Perform L-layer convolution operation on the time series data to obtain the deep feature map e L, e L The size of T*C L *W L *H L :

[0146] e l = ε l (e l-1 ) for l e [1, L]

[0147] where ε l () represents a convolution operation on the features of the l-1 layer, and when l = L, a convolution feature with a depth of L is obtained.

[0148] Step 302: The L-TAE multi-head self-attention module performs time series feature encoding on the feature map e L obtained in step 301, to obtain G groups of encoded time series feature weights where the size of each group of weights is T*W L *H L :

[0149]

[0150] Step 303: In order to use the time series feature weights on different depth encoders, the time series feature weights calculated in step 302 are resampled to have the same length and width as the shallow layer features:

[0151] a l,g = resize l (a L,g )

[0152] where resize l () represents resampling the feature in the spatial dimension to W l *H l .

[0153] Step 304: Apply the time series feature weights a l,g to the convolution feature map e l obtained by each layer of the encoder, and the application process is as follows:

[0154] Uniformly group each layer of features to obtain G groups of feature vectors with a size of T*C l / G*W l *H l : Multiply each feature vector by the time weight a l,g , and sum in the time dimension, and then perform a convolution operation on the result to obtain the time series encoded feature map:

[0155]

[0156] The time series data used in this part is collected in monthly units, and the image data of the target area for one year is collected. Due to the revisit period of Sentinel-2 and the cloud cover caused by the climate of the target area, it is generally impossible to collect complete 12-time image data of the target area in the actual implementation process, so only the image data with 12 time series is used for training, and the areas with incomplete time series are not involved in the model training.

[0157] During training, the parameters are set as follows:

[0158] The training epoch is 50, the batch size is set to 8, the initial learning rate is 0.001, Adam is used as the optimizer, and the weight decay is 1e -5 .

[0159] Additional training is performed on the model for difficult samples. Due to the existence of two or more products in some cultivated land samples, the recognition accuracy of the cultivated land product is low, and the recognition accuracy of another product is high. After actual image comparison, it is found that the first two products are wrong and the other product is correct. Therefore, during the preparation of weakly supervised labels, due to the inconsistency screening, the number of pixels of this type of cultivated land is very small, which causes the final model to not learn the characteristics of this type of cultivated land well, so some difficult sample data of the product are manually selected, about 30 samples of the cultivated land of this type are selected from the product with good prediction to continue learning, and the model is trained again. During training, since the product is directly used as the label, there is no need to add a weakly supervised label, and the training is set as follows: the epoch is set to 20, the batch size is set to 8, the initial learning rate is 0.001, Adam is used as the optimizer, and the weight decay is 1e -5

[0160] The model is used for cultivated land coverage mapping within a proper range. In the mapping link, the image data with incomplete time series can be used for cultivated land coverage prediction, and the 256x256 size samples obtained by sliding cropping are used as the unit for prediction. During the mapping process, the part of the three products with high consistency of cultivated land class, i.e., the part with a grid value of 1 in the weakly supervised data set, is retained, and the result of the weakly supervised model prediction is used to modify the area with inconsistent products, and finally the complete cultivated land coverage mapping result is obtained.

[0161] To illustrate the effectiveness of the present application, Hunan Province is taken as the research area, Sentinel-2 multispectral data is taken as the image data source, DynamicWorld, WorldCover of the European Space Agency and land cover products produced by the Esri company are taken as the cultivated land product data source, and the basic information of the three products is as follows:

[0162] Table 1: Product basic information

[0163] Product resolution Data source Model type Product time range Dynamic World 10m Sentinel 2 Self-built network 2017- World Cover 10m Sentinel 2 Unet 2017-2021 Esri 10m Sentinel 2 + DSM, etc. Random forest 2020-2021

[0164] The land use map obtained by Hunan Province 2 is used as the product precision verification standard, and the precision verification experiment is carried out on the data of 2020. The experimental process and results are as follows:

[0165] Hunan Province is uniformly divided into 121 block data units according to 0.5° latitude and longitude interval, and the data in each unit range is downloaded and processed. First, the remote sensing time series data missing degree is selected in the above unit range, and the data amount of the data set and the data set evenly covering the two requirements in the study area are considered. About 30 units with 0-2 months of missing time series are selected in the whole study area, and the corresponding time series image data and cultivated land product data are downloaded in each region. According to the weak supervision time series data set production method in step 1, 15488 data samples of 256x256 size are obtained after grid clipping. 3 / 5 of the data are used for weak supervision model learning, and the remaining part is used for model precision verification.

[0166] In order to prove that the method can produce higher precision cultivated land products based on existing products, the data of Hunan 2 is used as the true cultivated land verification standard, and the precision of multiple product sources and the final results obtained by the method are compared. The following four indicators are used for precision verification:

[0167] Precision (Precision), used to measure the result of the correct rate, that is, how many pixels predicted as cultivated land are correct:

[0168]

[0169] Recall (Recall), used to measure the result of the correct rate, that is, the ratio of all pixels that are truly cultivated land:

[0170]

[0171] F1 value is the weighted average of precision and recall, which comprehensively compares the two indicators:

[0172]

[0173] mIoU value, intersection over union:

[0174]

[0175] Among them, TP represents the number of pixels whose prediction result is correct and whose true result is positive class; FN represents the number of pixels whose prediction result is wrong and whose true result is positive class; FP represents the number of pixels whose prediction result is wrong and whose true result is negative class.

[0176] The final comparison accuracy is as follows:

[0177] Table 2. Comparison of different products and the mapping accuracy of the method

[0178]

[0179] The method proposed in the present application is higher than the existing product accuracy in three indicators except for the recall rate. In fact, since the method in the present application actually has no label cost, the time series image of Sentinel-2 can be obtained for free, and the utae model used in the present application can be trained at any time length, so theoretically, the model performance can be continuously enhanced by continuously adding weakly supervised labels on the basis of the current model.

[0180] The beneficial effects of the present application are as follows:

[0181] (1) The present application integrates the idea of weak supervision into the cultivated land product fusion mapping framework, solving the problem of a large number of effective pixels lost due to inconsistency between products in the old product fusion framework. And in view of the existence of a large number of existing free products and free remote sensing images, the present method can theoretically be promoted in global cultivated land mapping.

[0182] (2) For the problem of cultivated land thematic mapping, considering the seasonal change of cultivated land over time, the present application proposes a cultivated land product fusion mapping framework based on time series image data. In the framework, Sentinel-2 time series data is used as the image data source, and the utae model modified by the weak supervision module is used to encode the time series features of the cultivated land image, so that the mapping result is more in line with the time series characteristics of the cultivated land, effectively improving the accuracy of the cultivated land type in the previous land cover product.

[0183] The word "preferred" is used herein as a term of art to denote using, by way of example, illustration, or exemplification. Any aspect or design described herein as "preferred" is not necessarily to be construed as being more advantageous than other aspects or designs. The use of the word "preferred" is intended to present a choice or a preference only. The term "or" as used in this application is intended to mean an inclusive "or" rather than an exclusive "or". That is, unless specified otherwise, or as is clear from the context, the phrase "X employs A or B" is intended to mean any of the natural inclusive permutations. That is, the phrase "X employs A or B" is satisfied by any of the following instances: X employs A; X employs B; or X employs both A and B.

[0184] Moreover, although the present disclosure has been illustrated and described with respect to one or more implementations, equivalent alterations and modifications will occur to others skilled in the art based on the foregoing description and accompanying drawings. The present disclosure includes all such modifications and alterations and is limited only by the scope of the following claims. In particular regard to the various functions performed by the above described components (e.g., elements, engines, modules, etc.), the terms used to describe such components are intended to correspond, where appropriate, to any component which performs the specified function (e.g., that is functionally equivalent), even though not structurally equivalent to the disclosed structure which performs the function in the herein illustrated exemplary implementations of the present disclosure. In addition, while a particular feature of the disclosure can have been disclosed with respect to only one of several implementations, such feature can be combined with one or other features of the other implementations as can be desired and advantageous for any given or particular application. Furthermore, to the extent that the terms "including", "includes", "having", "has", "contain", "contains", or variants thereof to be afforded similar meanings in the context of describing common ownership in the background portion of the detailed description are used, such terms are intended to be inclusive in a manner similar to the term "comprising" as comparable terms under the doctrine of equivalents in the field of patent law.

[0185] The various functional units in the embodiments of the present application can be integrated in one processing module, or each unit can exist physically, or a plurality of or more units can be integrated in one module. The integrated module can be realized in the form of hardware, or in the form of a software functional module. When the integrated module is realized in the form of a software functional module and sold or used as an independent product, it can also be stored in a computer readable storage medium. The storage medium mentioned above can be a read-only memory, a magnetic disk or an optical disk, etc. The above-mentioned devices or systems can execute the storage method in the corresponding method embodiments.

[0186] In summary, the above-mentioned embodiments are one embodiment of the present application, but the embodiments of the present application are not limited by the above-mentioned embodiments, and any changes, modifications, substitutions, combinations and simplifications made without departing from the spirit and principle of the present application are equivalent replacement methods and are included in the protection scope of the present application.

Claims

1. A remote sensing farmland product fusion mapping method driven by weakly supervised temporal features, characterized in that, Includes the following steps: Construction of weakly supervised time series datasets; Construct a UTAe temporal network model modified with weakly supervised loss; Weakly supervised model training is used to obtain an initial model, followed by supplementary training for difficult samples. Finally, the trained model is used to create an appropriate area of ​​farmland cover mapping. The construction of the UTAe temporal network model modified by the weakly supervised loss includes: Step 201: Obtain the output of the last layer of the UTAE decoding layer, that is, the decoded image feature vector of the layer before the softmax layer of the network. At this time, the size of the feature vector is ,in For spatial dimension, In the channel dimension; in the spatial dimension, feature pixels are randomly selected at a rate of 1% to obtain a feature set of 0.01*W*H. , The size of each feature is ; Step 202: Calculate spatial dimensional correlation: For the feature set S obtained in step 201, perform the following operation on each element: Take this feature p i If an element is located on an image boundary, it is only 3-5 pixels within a spatial dimension of 8 adjacent pixels. Calculate the Euclidean distance between the element and the features of its spatially adjacent pixels, and select the two closest and farthest results as the two feature vectors closest and farthest from the pixel. Step 203: Calculate feature dimension correlation: For the feature set S obtained in step 201, perform the following operation on each element: Take this feature p i Multiple feature vectors in the channel dimension Construct a new set of channel features The size of each vector is Calculate in sequence The minimum cosine similarity among all elements is taken as the minimum similarity among the channel features. ; Step 204: Construct the weakly supervised loss: For each element in the feature set S, calculate the minimum and maximum spatial similarity distances of each element obtained in steps 202 and 203. and The minimum feature similarity calculated on the channel dimension Summing the three eigenvalues ​​yields the weakly supervised loss: in These are the weight parameters for each feature similarity, used to adjust the contribution of each similarity to the final loss; Step 205: Construct the complete loss of the final model: Calculate only the portion with label values ​​of 0 and 1, ignoring pixels with a raster value of 2 on the label. Cross-entropy calculation formula: in For the raster value corresponding to the label, The probability of the network prediction result; The final loss function formula is as follows: in Used to weight the weakly supervised loss, which is calculated based on the model's feature vector after the decoding layer.

2. The remote sensing farmland product fusion mapping method driven by weakly supervised temporal features according to claim 1, characterized in that, The process of creating weak supervision labels in the weak supervision time series dataset is as follows: For the given three data source products By performing tag mapping, all land types in the product except for cultivated land are converted into the same background pixels, while cultivated land pixels are retained, resulting in a transformed cultivated land-themed product. ; Add the three farmland thematic products together to obtain a single raster product with raster values ​​ranging from 0 to 3. For grid Perform label mapping, mapping pixels with a raster value of 1 to 2, and pixels with a raster value of 3 to 1, to obtain the final weakly supervised labels. In the above weakly supervised labels, the part with a grid value of 2 does not participate in the model training, the part with a grid value of 1 is the cultivated land pixel, and 0 is the background pixel; The formulaic description of this step is as follows: in Represents the original product. This is the first step of the label mapping rules. Represents grid addition. This is the label mapping rule for the second step.

3. The remote sensing farmland product fusion mapping method driven by weakly supervised temporal features according to claim 1, characterized in that, Discrete weakly supervised labels and Sentinel-2 time series data are used to supervise the learning of the modified UTAE time series model; The basic UTAE network adds an L-TAE module to the U-Net network for temporal feature encoding. L-TAE is a simplified multi-head self-attention mechanism network that directly processes the deep abstract features of the image after convolutional encoding, thereby transforming the basic U-Net model into a temporal UTAE model that can process time series data.

4. The remote sensing farmland product fusion mapping method driven by weakly supervised temporal features according to claim 3, characterized in that, The temporal feature encoding steps are as follows: Step 301: Perform L layers of convolution on the time-series data to obtain the depth feature map of the original image. , The size is : in Indicates to Convolution operations are performed on the features of the layer when This yields convolutional features of depth L. Step 302: The L-TAE multi-head self-attention module processes the feature map obtained in step 301. Perform temporal feature encoding to obtain G sets of encoded temporal feature weights. The size of each weight group is 1. : Step 303: In order to use temporal feature weights on encoders at different depths, the temporal feature weights calculated in step 302 are resampled to have the same width and height as the shallow layer features: in This indicates that features are resampled in the spatial dimension. , These are the weights of the sampled temporal features; Step 304: Weight the time series features Convolutional feature maps obtained by applying the algorithm to each encoder layer The process is as follows: Each layer of features is uniformly grouped to obtain G groups with a size of [missing value]. eigenvectors Using time weights For eigenvectors Perform term-by-term multiplication and summation over time, then perform a convolution operation on the result to obtain the temporally encoded feature map: 。 5. The remote sensing farmland product fusion mapping method driven by weakly supervised temporal features according to claim 4, characterized in that, Farmland cover mapping is performed using a model. In the mapping process, farmland cover prediction is performed using incomplete time-series image data. During mapping, 256×256 samples obtained by sliding cropping are used as the unit for prediction. The parts with high consistency of farmland land types among the three products are retained, that is, the parts with raster values ​​of 1 in the weakly supervised dataset. The prediction results of the weakly supervised model are used to correct the inconsistent areas of the products. Finally, the complete farmland cover mapping result is obtained by stitching together the data.

Images