一种基于时间序列遥感影像自动提取珊瑚礁的方法

专利名称：一种基于时间序列遥感影像自动提取珊瑚礁的方法

专利类型：发明专利

专利申请号：CN202110516844.1

专利申请（专利权）人：南京大学
权利人地址：江苏省南京市栖霞区仙林大道163号

专利发明（设计）人：程亮,庄启智,李嘉芮,陈辉,宋艳若,李满春,楚森森,东野升鹍,李泽明,刘磊

专利摘要：本发明涉及一种基于时间序列遥感影像自动提取珊瑚礁的方法，包括以下步骤：第一步、遥感影像并行预处理——对于影像进行大气校正；第二步、遥感影像自动筛选——从空间重叠、日期唯一、云量、图像熵四个方面实现影像的自动筛选；第三步、时间序列构建——构造影像MNDWI的时间序列；第四步、珊瑚礁自动提取——构建珊瑚礁时间序列的特征曲线，计算像素级时间序列与所述特征曲线之间的DTW值，使用二分法确定DTW阈值并提取珊瑚礁。本发明解决了珊瑚礁影像中存在多种噪声的问题，实现了遥感影像自动筛选，提出了一个可靠的基于时间序列遥感影像自动提取珊瑚礁的方法，并为基于其他卫星传感器自动提取珊瑚礁范围提供了流程思路。

主权利要求：
1.一种基于时间序列遥感影像自动提取珊瑚礁的方法，包括以下步骤：步骤1、遥感影像预处理——对遥感影像进行大气校正；
步骤2、遥感影像自动筛选——从空间重叠、日期唯一、云量、图像熵四个方面实现影像的自动筛选；其中，
1)、筛选去除空间重叠程度小的影像；
2)、通过遍历方法去除具有相同轨道号、相同影像获取日期和不同行列号的重复影像；
3)、筛选去除云覆盖率大的影像；
4)、筛选去除总图像熵大的影像；
步骤3、时间序列构建——基于MNDWI指数构建每个像素的时间序列，得到像素级时间序列；
MNDWI指数的计算方法如下：
其中，Green为绿光波段的灰度值，SWIR为短波红外波段的灰度值；
步骤4、珊瑚礁自动提取——选择任意一幅成像质量好的影像，获取珊瑚礁真实面积，从该幅影像中随机选取若干像元，通过目视判别将其进行分类，分为珊瑚礁或海洋，分别将每幅影像中与珊瑚礁坐标点对应的像元的MNDWI指数进行平均，使得每幅影像得到一个关于珊瑚礁的MNDWI指数，进而获得关于珊瑚礁的基于MNDWI指数的特征曲线，作为珊瑚礁时间序列的真值；计算像素级时间序列与所述特征曲线之间的DTW值，使用二分法确定DTW阈值使得提取的珊瑚礁面积与珊瑚礁真实面积之比达到预设阈值，将DTW值小于DTW阈值的像元划分为珊瑚礁并进行提取。
2.根据权利要求1所述的基于时间序列遥感影像自动提取珊瑚礁的方法，其特征在于：步骤2中，空间重叠程度利用裁剪影像的有效范围百分比来表示，其比例计算如下：其中，areaeffective代表有效面积，不包括黑色边缘和无数据的区域，areaall代表实验区域的面积，pspace表示空间重叠百分比，基于空间重叠进行自动筛选的阈值取95％。
3.根据权利要求1所述的基于时间序列遥感影像自动提取珊瑚礁的方法，其特征在于：步骤2中，云覆盖率计算如下：
其中，areacloud代表检测出云的面积，areaall代表实验区域的面积，pcloud表示云覆盖率；
基于云量进行自动筛选的阈值取5％。
4.根据权利要求1所述的基于时间序列遥感影像自动提取珊瑚礁的方法，其特征在于：步骤2中，分别计算蓝光波段、绿光波段和近红外波段的图像熵，图像熵H计算公式如下：2
pij＝f(i，j)/N
其中i是像素的灰度值，j是邻域的灰度值，f(i，j)是特征二元组(i，j)的频率，N是图像的尺度，pij反映周围像素灰度分布的综合特征；
计算影像的总图像熵Hsum
Hsum＝HBlue+HGreen+HNIRHBlue为蓝光波段图像熵，HGreen为绿光波段图像熵，HNIR为近红外波段的图像熵，总图像熵的阈值取5。
5.根据权利要求1所述的基于时间序列遥感影像自动提取珊瑚礁的方法，其特征在于：步骤4中，所述珊瑚礁真实面积通过目视选取或先验数据获得。
6.根据权利要求1所述的基于时间序列遥感影像自动提取珊瑚礁的方法，其特征在于：步骤4中，所述DTW值替换为欧氏距离、模式距离或形状距离。 说明书 : 一种基于时间序列遥感影像自动提取珊瑚礁的方法技术领域[0001] 本发明涉及一种基于时间序列遥感影像自动提取珊瑚礁的方法，特别是涉及一种可以实现遥感影像自动筛选、基于时空相似性自动提取珊瑚礁的方法。背景技术[0002] 珊瑚礁是地球上最具多样化和经济价值的生态系统之一。珊瑚礁仅覆盖地球表面的一小部分(不到百分之一)以及不到百分之二的海床。然而，有四分之一的海洋物种依靠珊瑚礁来提供食物和庇护所，因此珊瑚礁通常被称为海洋雨林(Mulhall，2009；Bellwood，2004)。珊瑚礁还可以提供一些重要的生态系统服务，例如渔业资源、沿海保护和旅游业。受自然变化(全球变暖、海平面上升、海洋酸化等)和人类活动(过度捕捞、环境污染、港口等基础设施等)的影响，珊瑚礁在全球范围内退化(Jackson，2001；Pandolfi，2003；Hoegh‑Guldberg，2007；De'ath，2012；Glynn；2015；McClenachan，2017；Barkley，2015)。对于珊瑚礁的监测与管理，基线信息十分重要。因此，迫切需要生成公共领域的珊瑚礁地图集，以有效指导珊瑚礁保护(Goodman，2013；Hamylton，2017)。如果没有这些基本的基线信息，就很难在全球范围内对珊瑚礁危机作出反应。同时，珊瑚礁基线形成的深水掩膜对于提高珊瑚礁分类精度也有着重要作用(Zainal,1993；Reshitnyk，2014)。[0003] 为了更好地测量珊瑚礁基线范围，学者们最开始主要采用实地测量和多数据源融合的方式。NOAA对太平洋和加勒比海海域0至150米深度范围的43000平方公里海域进行了测绘，但是实地测量不能覆盖全部珊瑚礁，只能局限于珊瑚礁的某一部分，而且争议海域的实地测量存在困难(Monaco，2012)。为了提供全球珊瑚礁的位置和现状信息，联合国环境规划署世界保护监测中心(UNEP‑WCMC)利用遥感影像和近150年的底图汇编成一张世界珊瑚礁地图集，并在2018年生产了迄今为止最全面的全球温水珊瑚礁数据集，为将来开展更详细的统计工作提供了基础基线图(Spalding，2001；UNEP‑WCMC‑Global，2019；Andrefouet，2006)。多源数据融合的方式虽然可以绘制出大范围区域的珊瑚礁基线图，但是由于不同来源地图的制图方法不同，存在准确性不一致的问题。[0004] 成熟的卫星遥感技术为珊瑚礁大范围周期性的观测提供了机会。利用遥感影像对珊瑚礁进行测绘的方法可以归纳为手动数字化、基于像素、面向对象(OBIA)三种。一些学者基于Landsat7ETM+影像手动提取珊瑚礁，并计算了珊瑚礁的面积(Naseer，2004)。先前的研究还证实了可以通过融合激光、声呐、剖面、航空摄影数据手动绘制0‑35米水深内的珊瑚礁类别(Walker，2008)。但是手动数字化是主观的，且耗费大量时间和精力。为了在不牺牲精度的情况下使方法自动/半自动化，基于像素(最大似然、光谱指数、深度学习等)和面向对象(SVM、随机森林、决策树、KNN等)方法开始应用于珊瑚礁的分类与提取。利用最大似然分类方法(基于像素)对珊瑚礁高分辨率影像(worldview、QB、IKONOS等)进行分类，分类精度可以达到71‑75％(Reshitnyk，2014；Capolsini，2003；Rowlands，2008；Collin，2012；Naidu，2018；Gholoum，2019；Phinn，2012；Leon，2011；Ahmad，1994；Dong，2019)。在面向对象方法中，SVM方法被证明对珊瑚礁分类效果最好，精度最高可以达到78％(Gholoum，2019；Wahidin，2015；Baumstark，2016；Selgrath，2016；Kabiri，2018；Roelfsema，2018；Purkis，2019)。随着深度学习在遥感图像分析领域的应用，深度学习和集成学习的方法逐渐流行。通过集成SVM、KNN、随机森林等方法，结合高光谱、航空摄影和水深测量等多种数据源，珊瑚礁分类精度提高到88.5％(Zhang，2015)。通过CNN方法训练PlanetDove高分辨率卫星制作珊瑚礁样本，珊瑚礁提取精度可以达到93％(Li，2020)。[0005] 遥感影像的分辨率和影像质量对珊瑚礁基线提取至关重要。影像分辨率越高，珊瑚礁提取精度越高(Gholoum，2019；Andrefouet，2003)。相较于航空影像，卫星遥感影像具备范围大、重复性观测、成本低等优势(Andrefouet，2006；Purkis，2019；Mumby，1997；Yamano，2013)。海洋地区遥感影像常受云、浪花、太阳耀斑等噪声的影响，限制了珊瑚礁提取精度(Reshitnyk，2014；Polcyn，1976；Polcyn，1970)。一些学者结合了大气校正、云遮罩和低通滤波器进行图像降噪(Traganos，2018)。但是，这些预处理技术很难实现大规模的一致应用(Chen，2004；Hedley，2016)。Sentinel‑2影像的云提取有很多成熟和开源的算法，比如S2cloudless(Zupanc，2019)、MAJA(Hagolle，2017)、Sen2Cor(Main‑Knorn，2017)、Fmask(Qiu，2019；Zhu，2015；Qiu，2017；Frantz，2018)等。相较之下，Fmask4.0具备较高的精度，精度可以达到90％(Sanchez，2020；Baetens，2019)。光学多波段影像的时间序列可提供长期的特征信息，这为大规模绘制珊瑚礁地图提供了机会(Wang，2019)。保留多时相影像中的无云像素可以构建一幅无云、无浪花的珊瑚礁影像(Chu，2019)。[0006] 显然，必须发展一种提取珊瑚礁的方法，解决珊瑚礁影像中存在多种噪声的问题，实现遥感影像自动筛选，并基于时间序列遥感影像对于珊瑚礁进行自动提取。发明内容[0007] 本发明要解决的技术问题是：克服现有技术的上述不足，提出一种基于时间序列遥感影像自动提取珊瑚礁的方法。相对于传统的方法而言，本发明充分考虑空间、时间、云等因素，实现了遥感影像自动筛选，解决了珊瑚礁影像中存在多种噪声的问题，提出了一种基于MNDWI指数和像素时空相似性分析自动提取珊瑚礁的方法，为珊瑚礁的提取提供了一个科学的实施方法。[0008] 为解决上述技术问题，本发明使用二分法确定DTW阈值，使得提取的珊瑚礁面积与珊瑚礁真实面积之比达到预设阈值，将DTW值小于DTW阈值的像元划分为珊瑚礁并进行提取。[0009] 本发明提出的珊瑚礁提取方法有效效益如下：[0010] (1)构造时间序列需要大量影像，本发明通过从空间重叠、日期唯一、云量、图像熵四个方面实现了影像的自动筛选，减少了人工选择影像的工作量；[0011] (2)本发明基于时空相似性分析，提出的框架能够以较高的准确性、完整性和正确性提取珊瑚礁，相较于从单一场景图像中提取珊瑚礁的其他方法，该框架可以有效地消除噪声并完全提取出珊瑚礁；[0012] (3)本发明运用开源数据构建时间序列，未来可以对大范围珊瑚礁进行提取，便于珊瑚礁的监测与管理。[0013] 本发明解决了珊瑚礁影像中存在多种噪声的问题，实现了遥感影像的自动筛选，提出了一个可靠的基于时间序列遥感影像自动提取珊瑚礁的方法，并为基于其他卫星传感器自动提取珊瑚礁范围提供了流程思路。附图说明[0014] 下面结合附图对本发明作进一步的说明。[0015] 图1为本发明实例总体流程图。[0016] 图2为本发明实例研究区概况图。[0017] 图3为本发明实例遥感影像自动筛选示意图。[0018] 图4为本发明实例珊瑚礁真值时间序列构建图。[0019] 图5为本发明实例珊瑚礁提取结果图。[0020] 图6为本发明实例不同DTW阈值的精度评价图。具体实施方式[0021] 下面根据附图详细阐述本发明，使本发明的技术路线和操作步骤更加清晰。[0022] 安达礁位于南沙群岛郑和群礁的东部，是郑和群礁中最大的岛礁，中心经纬度为北纬10°20′，东经114°42′。如图2所示，安达礁位于太平岛以东36公里(20海里)，南薰礁东北约54公里(29海里)，面积约30平方公里。[0023] 考虑到安达礁面积较大且地形所处水深不同，选择安达礁作为研究区域，代表南沙群岛的大多数珊瑚礁。安达礁是一个未封闭的环礁，尖端朝向东北，钳状的部分朝向西南。安达礁的东北部相对较浅，最浅的部分靠近海平面，西南较深，水深超过30米。安达礁前缘是一个剑形的山脊，延伸到海中，礁脊的两侧陡峭而狭窄，水深从9.1米下降到91.4米，这对安全航行构成了威胁。[0024] 本实例以该研究区为例对一种基于时间序列遥感影像自动提取珊瑚礁的方法进行说明，如图1流程图所示，具体包括以下步骤：[0025] 步骤1、遥感影像并行预处理——在Linux集群建立用于大气校正的Docker环境，利用Sen2Cor(参见论文Main‑Knorn,M.,etal.,Sen2CorforSentinel‑2.ImageandSignalProcessingforRemoteSensing,2017.10427.)对遥感影像进行大气校正。[0026] 安达礁地区2015年11月26日至2020年12月24日Sentinel‑2卫星影像共598景。将598个Sentinel‑2影像分为20组，在Linux集群(128核心处理器，1856GB内存，298T存储空间)中构建用于大气校正的docker环境。在Pipeline技术的支持下，并行处理大约需要6.5个小时。由于早期和当前Sentinel‑2影像的数据组织方法不同，因此需要使用不同版本的Sen2Cor进行大气校正。如果使用14.2产品标准生成Sentinel‑2Level‑1C产品数据，大气校正需使用Sen2Cor_v2.5.5版本，其他产品的大气校正使用最新版本Sen2Cor_v2.8。[0027] 步骤2、遥感影像自动筛选——本发明主要从空间重叠、日期唯一、云量、图像熵四个方面实现影像的自动筛选，如图3所示。[0028] 基于空间重叠的自动筛选即筛选去除空间重叠程度较小的影像。空间重叠程度利用裁剪影像的有效范围百分比来表示，其比例计算如下：[0029][0030] 其中，areaeffective代表有效面积，不包括黑色边缘和无数据的区域，areaall代表实验区域的面积，pspace表示空间重叠百分比(0≤pspace≤1)。[0031] 基于日期唯一的自动筛选即通过遍历方法去除具有相同轨道号、相同影像获取日期和不同行列号的重复影像。[0032] 基于云量的自动筛选即通过Fmask(参见论文Qiu,S.,Z.Zhu,andB.B.He,Fmask4.0:ImprovedcloudandcloudshadowdetectioninLandsats4‑8andSentinel‑2imagery.RemoteSensingofEnvironment,2019.231.)对于云覆盖率进行评估，并筛选去除云覆盖率过大的影像，云覆盖率计算如下：[0033][0034] 其中，areacloud代表检测出云的面积，areaall代表实验区域的面积，pcloud表示云覆盖率(0≤pcloud≤1)。[0035] 基于图像熵的自动筛选即筛选去除图像熵较大的影像，可以选择图像的邻域灰度平均值作为灰度分布的空间特征来计算二维熵(参见论文Corominas‑Murtra,B.,J.Fortuny,andR.V.Sole,Towardsamathematical theoryofmeaningfulcommunication.ScientificReports,2014.4.)。[0036] 首先分别计算蓝光波段、绿光波段和近红外波段的图像熵，图像熵H计算公式如下：[0037][0038] pij＝f(i,j)/N2[0039] 其中i是像素的灰度值(0≤i≤255)，j是邻域的灰度值(0≤j≤255)，f(i,j)是特征二元组(i,j)的频率，N是图像的尺度，pij可以反映周围像素灰度分布的综合特征,H为图像的图像熵，反映了图像的复杂度，蓝光波段和绿光波段在水中具有很强的透射能力，而近红外波段在云中具有很强的穿透能力，故选择这些波段来计算图像的总图像熵Hsum：[0040] Hsum＝HBlue+HGreen+HNIR[0041] HBlue为蓝光波段图像熵，HGreen为绿光波段图像熵，HNIR为近红外波段的图像熵。如果Hsum小，则表明图像中的信息相对稳定，并且气泡云、云阴影和波浪很少，否则应删除图像。[0042] 将经过图像预处理的598景Sentinel‑2影像通过以下约束条件进行筛选：[0043] (1)空间重叠百分比大于95％；[0044] (2)影像获取日期唯一；[0045] (3)云量小于5％；[0046] (4)在安达礁范围内的Sentinel‑2影像总图像熵不能超过5。[0047] 在这种情况下，安达礁研究区域的Sentinel‑2影像被筛选至127景。[0048] 由于Sentinel‑2卫星从2015年底才开始获取中国南海的影像，因此进一步分析2016年至2020年每个季节的影像数量。可以发现第一季度保留的影像最多，而第二、第三季度保留的影像较少，该统计规则也与南海热带季风气候的雨季相一致，即冬春季节少雨，夏秋季节多雨。通过对统计信息的分析和验证，可以证明所提出自动筛选影像的可行性和稳定性。[0049] 步骤3、时间序列构建——为了突出珊瑚礁在时间序列构建中的特征，选择安达礁2019年3月30日的Sentinel‑2影像，分析其绿光波段、近红外波段和短波红外波段的特征，并进行定性比较。经分析得出，MNDWI可以更好地显示出珊瑚礁的特征。将获得的MNDWI结果按时间顺序排序，在图像严格对齐之后构造每个像素的时间序列；其中，MNDWI的计算方法如下：[0050][0051] 步骤4、珊瑚礁提取——选择任意一幅成像质量较好的影像，获取珊瑚礁真实面积(通过目视选取或先验数据获得)，随机生成100个像素点，通过目视判别将其划分为珊瑚礁与海洋，分别将每幅影像中与珊瑚礁坐标点对应的像元的MNDWI指数进行平均，使得每幅影像得到一个关于珊瑚礁的MNDWI指数，进而获得关于珊瑚礁的基于MNDWI指数的特征曲线，作为珊瑚礁时间序列的真值。通过DTW(参见论文Chen,T.S.,T.Y.Lin,andY.W.P.Hong,GaitPhaseSegmentationUsingWeightedDynamicTimeWarpingandK‑NearestNeighborsGraphEmbedding.2020IeeeInternationalConferenceonAcoustics,Speech,andSignalProcessing,2020:p.1180‑1184.)分析像素级时间序列的相似性(即，计算像素级时间序列与所述特征曲线之间的DTW值)，使用二分法确定DTW阈值使得提取的珊瑚礁面积与珊瑚礁真实面积之比达到预设阈值，将DTW值小于DTW阈值的像元划分为珊瑚礁并进行提取。[0052] 在对安达礁进行2015年至2020年时间序列的构建之后，有必要弄清楚属于珊瑚礁的像素时间序列的曲线特征。根据2020年9月20日的Sentinel‑2影像在实验区域中随机选择100个采样点，将其手动分为珊瑚礁或海洋，有42个采样点属于珊瑚礁，58个采样点位于海洋。分别将每幅影像中与珊瑚礁坐标点对应的像元的MNDWI指数进行平均，使得每幅影像得到一个关于珊瑚礁的MNDWI指数，使用S‑G滤波降噪方法(参见论文Chen,J.,etal.,Asimplemethodforreconstructingahigh‑qualityNDVItime‑seriesdatasetbasedontheSavitzky‑Golayfilter.RemoteSensingofEnvironment,2004.91(3‑4):p.332‑344.)进行平滑，进而获得关于珊瑚礁的基于MNDWI指数的特征曲线，作为珊瑚礁时间序列的真值，如图4所示。[0053] 利用珊瑚礁时间序列的真值来测量Sentinel‑2影像构建的像素级时间序列的相似度。在时空相似性测量结果中，最大DTW为2027.4111，最小为0.8572。经过二分法自动选择阈值11轮之后，最终确定DTW阈值为25.5954，相应的安达礁提取结果如图5所示。当阈值为10和40时，安达礁提取结果如图6所示，前者仅提取出了部分珊瑚礁，后者存在错误提取的部分，故当阈值DTW为25.5954时，安达礁具有最佳提取结果。[0054] 下面为验证本发明方法的精度及可靠性，继续以该实例进行说明。提取结果分为珊瑚礁和非珊瑚礁，手动将无云的珊瑚礁影像以及珊瑚礁提取结果数字化，并创建混淆矩阵，计算总体准确性、精确性、召回率和F1分数进行精度评价。当阈值DTW为25.5954时，安达礁具有最佳提取结果，分别为98.59％、98.04％、96.05％和0.9703。[0055] 为了验证基于时间序列遥感影像的珊瑚礁提取的可靠性和稳定性，选择大津算法(参见论文Otsu,N.,ThresholdSelectionMethodfromGray‑LevelHistograms.IeeeTransactionsonSystemsManandCybernetics,1979.9(1):p.62‑66.)、人工分割、最大似然法(参见论文Reshitnyk,L.,etal.,EvaluationofWorldView‑2andacousticremotesensingformappingbenthichabitatsintemperatecoastalPacificwaters.RemoteSensingofEnvironment,2014.153:p.7‑23.)、支持向量机(参见论文Zhang,C.,Applyingdatafusiontechniquesforbenthichabitatmappingandmonitoringinacoralreefecosystem.IsprsJournalofPhotogrammetryandRemoteSensing,2015.104:p.213‑223.)和MaskR‑CNN(参见论文He,K.M.,etal.,MaskR‑CNN.2017IeeeInternationalConferenceonComputerVision(Iccv),2017:p.2980‑2988.)在相同实验环境下从无云Sentinel‑2图像(2020‑09‑20)中提取安达礁石，不同方法的总体准确性、精确性、召回率和F1分数如下表所示。[0056] 表1不同方法精度评价[0057][0058] 本发明安达礁提取准确性比大津算法高14.81％，比人工分割高5.53％，比最大似然法高7.52％，比支持向量机高1.65％，比MaskR‑CNN高％1.95。本发明安达礁提取的F1分数分别比大津算法、人工分割、最大似然法、支持向量机和MaskR‑CNN高出0.4309、0.1262、0.1973、0.0354和0.0385。人工分割的精度最高，比本发明的分类精度高0.93个百分点，而召回率比本发明的分类精度低22.46个百分点。可见本发明精度较高，可靠性很强。[0059] 由于Landsat8影像产品经过了几何校正，因此本发明提出的框架也具有一定的适用性。为了验证该方法的可扩展性，根据相同的实验过程，使用Landsat8影像提取了安达礁。[0060] 选择2015年1月30日至2020年12月29日的85张Landsat8影像。经过大气校正后，根据以下约束条件选择了52幅图像：[0061] (1)空间重叠百分比大于95％；[0062] (2)影像获取日期唯一；[0063] (3)云量小于5％；[0064] (4)在安达礁范围内的Sentinel‑2影像总图像熵不能超过5。[0065] 构建珊瑚礁时间序列的特征曲线，计算像素级时间序列与所述特征曲线之间的DTW值，使用二分法确定DTW阈值并提取珊瑚礁。在时空相似性测量结果中，最大DTW为3485.5302，最小为0.0051。经过二分法13轮自动阈值选择后，最终DTW为5.5395，此时安达礁具有最佳提取结果，总体准确性、精确性、召回率和F1分数分别为98.59％、96.71％、97.32％和0.9701，精度较高，可见本发明的可扩展性较好。[0066] 本发明基于时间序列遥感影像自动提取珊瑚礁的方法不局限于上述实施例所述的具体技术方案，凡采用等同替换形成的技术方案均为本发明要求的保护范围。

专利地区：江苏

专利申请日期：2021-05-12

专利公开日期：2024-07-26

专利公告号：CN113128523B