弱监督时序特征驱动的遥感耕地产品融合制图方法

本发明属于遥感耕地产品的制图，尤其涉及弱监督时序特征驱动的遥感耕地产品融合制图方法。

背景技术：

1、耕地作为人类生存和社会发展的基本资源，其面积和分布受到我国政府的密切关注，要严守18亿亩耕地红线，及时准确的耕地制图至关重要。随着遥感技术的快速发展，大批不同时空分辨率，光谱分辨率的卫星传感器的升空为更快速更高精度的地表资源监测提供了机会。目前，有两类产品能够为耕地监测服务，耕地专题制图产品与全球土地覆盖地图。

2、对于耕地专题图制作，各国开展了一系列相关项目，常见的耕地制图产品有美国的cdl(cropland data layer)作物分类图，加拿大的afcc产品以及欧洲的sen-agri，由于这些地图的制作往往需要制作足量的标签数据用于模型的学且，且需要大量的实地调查验证产品可靠性，因此产品覆盖范围只局限在上述项目实施地区，通常难以大范围的进行耕地专题地图的制作。

3、对于全球土地覆盖产品，不同国家以及组织发布了许多区域或全球范围的土地覆盖产品。例如，谷歌的10m近实时的全球地表覆盖项目：dynamic world，esri公司发布的2017-2021为期5年的地表覆盖产品，清华大学的from-glc等。这些高空间分辨率的地表覆盖产品的出现为全球尺度下的耕地监测提供了机会，但是由于这类产品首要保障产品全球覆盖且综合保障所有土地覆盖类型整体的精度，因此为保证遥感数据的全球覆盖，常使用多个影像拼接而成单一时相数据或时间间隔较长的数据进行分类。而耕地这一覆盖类型的光学特征随季节改变有明显的变化，相同地区不同时间的耕地影像在纹理，色彩等视觉特征也会发生很大改变，因此使用单一时相的遥感影像制图得到的土地覆盖结果中耕地的精度并不高。

4、由于耕地的分布和面积一般不会随时间有剧烈的变化，为了低成本的进行高精度的全球耕地覆盖产品制作和更新，许多研究者尝试将多个土地覆盖产品根据一定的标准融合来生成新的精度更高的耕地专题图。下面从产品融合的数据源和产品融合方法两方面来介绍当前的耕地产品融合技术现状与技术不足之处。

5、产品融合数据源：遥感影像的空间分辨率限制了土地覆盖产品的制图精度，早期土地覆盖产品以较低的空间分辨率为主，随着欧洲哨兵系列卫星升空，我国的高分、资源系列卫星相继开放共享数据，土地覆盖产品的制作迈进了中、高空间分辨率时代。随着机器学习方法和遥感云计算平台的发展，高空间分辨率，高频次的全球土地覆盖制图成为可能，当前常用的几个全球土地覆盖产品如下：

6、谷歌的dynamic world产品：谷歌公司在2022年发布，其空间分辨率为10m，覆盖全球范围，时间跨度为2015-2022年，使用哨兵2号的1c级产品作为影像，该产品使用深度学习模型作为分类算法，以24000个全球各地标注图像为训练基础，凭借谷歌自身的云平台高强度算力，做到了近实时的土地覆盖制图。

7、欧空局的world cover产品：esa在2021年发布的产品，wc产品发布了2020年与2021年两期，空间分辨率为10m，覆盖全球范围，wc使用哨兵2号和哨兵1号作为基础数据，综合使用了osm等开源地理数据作为辅助，用随机森林分类树算法在100×100米网格中手工标记的像素上进行训练预测得到。

8、esri的land cover产品：esri公司2021年发布的产品，空间分辨为10m，时间跨度为2017-2021，覆盖全球范围，使用哨兵2号作为模型输入数据，每期产品使用单时相数据。esri的产品使用深度学习方法制作，采用的训练数据集与谷歌的产品相同，但在分类结束后根据先验对分类结果做了微调。

9、然而，一方面，这些产品的制作过程没有考虑到耕地的时空异质性，另一方面，由于不同数据之间使用的土地覆盖分类体系不同，分类方法不同，验证方法不同，不同产品之间在协同融合时必然会出现较大不一致性，不同的产品直接叠加可能将不同产品的错误分类结果带到融合产品从而导致融合精度降低。

10、产品融合技术方法：产品融合方法需要克服现有土地覆盖产品自身缺陷以及不同产品之间的不协调、不兼容的限制，才能融合得到更高精度的融合结果。当前的产品融合方法可以从学习式和非学习式两类产品融合方法来分类介绍。

11、非学习式：非学习式方法试图结合统计的数据融合方法来解决现有产品之间兼容性、可比性和准确性不够的问题。按照融合方法这些研究可以分为基于地理加权回归(gwr)的和基于融合决策规则的。gwr方法其特征是回归参数随着空间距离的变化而变化，因此能够较好的反应空间位置因素带来的影响。linda see等人在非专利文献“building ahybrid land cover map with crowdsourcing and geographically weightedregression，isprs journal of photogrammetry and remote sensing,103:48-56，2015”通过gwr方法将glc2000,modis全球土地覆盖产品，globcover三个产品融合，生成了全球300m产品图。

12、融合决策规则包含有多种模型，贝叶斯理论，dempster-shafer证据理论。贝叶斯理论在已知先验概率分布的条件下可提高类别之间的可分性，a.pérez-hoyos等在非专利文献：integrating multiple land cover maps through amulti-criteria analysis toimprove agricultural monitoring in africa,international journal of appliedearth observation and geoinformation，2020”中使用贝叶斯定理制作了融合的全球土地覆盖地图，给出了每个像素分类的置信度值，减少了因数据协调性和空间采样等因素产生的误差。

13、学习式：学习式的产品融合方法试图收集到不同产品高可信度的样本，将其作为标签输入到分类器中，并通过额外的遥感影像重新生成一幅精度更高的局部或全球的地图。由于产品之间质量参差不齐，通常需要一个样本筛选机制来选取可靠样本，常用的方法有基于规则筛选，如只保留多年数据均没有发生变化的样本，只保留中心像素以减少边缘效应的影响，或使用边缘缩减算子以去除边界像素，liu等人在非专利文献“high-resolution multi-temporal mapping of global urban land using landsat imagesbased on the google earth engine platform，remote sensing ofenvironment,2018,209:227-239”通过gee(google earth engine)进行像素滤波和光谱滤波从modis地表覆盖产品中筛选标签，最后得到样本精度达到了99.2％。

14、然而，一方面，非学习式的数据融合方法虽然可以减少不同产品之间不协调、不兼容的问题，但是缺少遥感影像辅助，无法针对耕地专题图做专门优化，原始产品使用单时序得到的结果，在这里仍旧是单时序的结果。另一方面，为减少产品自身和产品之间的误差，现有的学习式融合框架需要抛弃不同产品之间的不一致像素，从而难以获取到完整且准确的深度学习所需标签，而深度学习的监督分类方法通常需要完整的准确标签作为模型指导才能训练得到可靠的分类模型，这也导致传统的深度学习方法不能直接用于学习式的产品融合方法，当前的学习方法均以机器学习的方法为主。最后，无论是非学习式还是学习式数据融合框架均没有考虑到耕地的时序性特，融合过程中仍使用单一时相的影像数据，最后精度没有得到进一步的提高。

技术实现思路

1、有鉴于此，本发明的目的是更高精度的融合当前土地覆盖产品以快速生产耕地覆盖产品。为此，一方面考虑将弱监督学习方法引入当前的产品融合框架中以改进学习式的产品融合方法，以在标签不完整的情况下学习得到可靠的模型用于耕地覆盖制图。另一方面，考虑使用长时序影像数据作为遥感数据源进行耕地时序特征建模，以utae网络，一种带时序特征编码结构的深度学习网络为基础模型，进行耕地时序特征建模以提升最终耕地覆盖制图精度。基于上述两点，本发明提出了一种基于弱监督学习时序特征提取的遥感耕地产品融合制图框架。

2、具体地，本发明提出的弱监督时序特征驱动的遥感耕地产品融合制图方法，包括以下步骤：

3、弱监督时序数据集构建；

4、构造由弱监督损失改造后的utae时序网络模型；

5、弱监督模型训练，得到初步的模型后再针对难样本进行补充训练；

6、最后使用训练好的模型进行适当范围的耕地覆盖制图。

7、进一步地，所述构造由弱监督损失改造后的utae时序网络模型，包括：

8、步骤201：获取utae解码层最后一层的输出，即网络softmax层前一层的经过解码后的影像特征向量vf，此时特征向量尺寸为c*w*h，其中w,h为空间维，c为通道维；在空间维度上，以1％的比例随机选取特征像素获得数量为0.01*w*h的特征集合s＝{p1,p2,p3...pn}，n＝0.01×w×h，其中每个特征的尺寸为c*1*1；

9、步骤202：计算空间维度相关性：对于步骤201中得到的特征集合s，对其中的每个元素做如下操作：

10、取该特征pi在空间维度上相邻的8个像素，若处于图像边界上则只有3-5个像素；计算该元素与其空间维度上相邻像元特征的欧式距离，并选取距离最近和距离最远的两个结果作为与该像素距离最近和最远的两个特征向量；

11、步骤203：计算特征维度相关性：对于步骤201中得到的特征集合s，对其中的每个元素做如下操作：

12、取该特征pi在通道维度上多个特征向量v'构成一个新的通道特征集合sv，其中每个向量的尺寸为1*1*1，依次计算sv中每个元素之间的余弦相似度，取其中的最小值作为通道特征上的最小相似度cmin；

13、步骤204：构造弱监督损失：对特征集合s中每个元素在步骤202与步骤203中对所做计算得到的每个元素在空间上最小和最大空间相似度距离dmin和dmax，在通道维度上计算得到的最小特征相似度cmin三个特征值做求和操作，得到弱监督损失：

14、

15、其中α,β,γ分别为各个特征相似度的权重参数，用于调整各个相似度对最终损失的贡献；

16、步骤205：构造最终模型的完整损失：只在标签值0和1的部分计算，忽略标签上栅格值为2的像元，交叉熵计算公式：

17、lossce＝σ-(yi·log(pi)+(1-yi)·log(1-pi))

18、

19、其中yi为标签对应栅格值，为网络预测结果概率；

20、最终的损失函数公式如下：

21、loss＝lossce+λ·lossweak

22、其中λ用于对弱监督损失加权，弱监督损失计算是基于模型的经过解码层之后的特征向量。

23、进一步地，使用离散的弱监督标签和哨兵2号时序数据来监督改造后的utae时序模型学习；

24、基础的utae网络在u-net网络的基础上添加l-tae模块进行时序特征编码，l-tae为一种简化的多头自注意力机制网络，直接处理图像经过卷积编码后的深度抽象特征，从而将基础的u-net模型改造为可处理时间序列的时序utae模型。

25、进一步地，所述时序特征编码步骤如下：

26、步骤301：对时序数据经过l层卷积操作，得到原始图像的深度特征图el，el的尺寸为t*cl*wl*hl：

27、el＝εl(el-1)for l∈[1,l]

28、其中εl()表示对l-1层的特征进行卷积操作，当l＝l时得到深度为l的卷积特征；

29、步骤302：l-tae多头自注意力模块对步骤301得到的特征图el进行时序特征编码，得到g组编码后的时序特征权重其中每组权重的尺寸均为t*wl*hl：

30、

31、步骤303：为了在不同深度的编码器上使用时序特征权重，从步骤302计算得到的时序特征权重重采样至长宽与浅层特征相同：

32、al,g＝resizel(al,g)

33、其中resizel()表示在空间维度将特征重采样至wl*hl，al,g为采样后的时序特征权重；

34、步骤304：将时序特征权重al,g作用于每一层编码器得到的卷积特征图el，作用过程如下：

35、将每一层特征均匀分组，得到g组尺寸为t*cl/g*wl*hl的特征向量使用时间权重al,g对特征向量做逐项乘法，并在时间维度上求和，对结果做一次卷积操作得到时序编码后的特征图：

36、

37、进一步地，使用模型进行耕地覆盖制图，在制图环节，使用时序残缺的影像数据进行耕地覆盖预测，制图时使用滑动裁剪得到的256×256大小的样本为单位预测；保留三个产品中耕地地类一致性较高的部分，即弱监督数据集中栅格值为1的部分，使用弱监督模型预测的结果来修改产品不一致的区域，最终拼接得到完整的耕地覆盖制图结果。

38、本发明的有益效果如下：

39、(1)本发明将弱监督的思想融入耕地产品融合制图框架中，解决了旧有产品融合框架中由于产品之间不一致性带来的大量有效像素丢失的问题。且鉴于存在大量已有免费产品和免费的遥感影像，因此本方法理论上可以在全球范围内进耕地制图推广。

40、(2)针对耕地专题制图问题，考虑到耕地随时间的季节性变化问题，本发明提出了以时序影像数据为基础的耕地产品融合制图框架，框架中采用哨兵2号时序数据为影像数据源，使用弱监督模块改造后的utae模型对耕地影像进行时序特征编码，让制图结果更符合耕地的时序性特征，有效的改良了以往土地覆盖产品中耕地类型的精度。