一种提取南黄海浒苔的自适应阈值遥感算法

王怡人，王胜强，2，喻樾，焦红波，孙德勇，2

(1.南京信息工程大学 海洋科学学院，南京 210044；2.江苏省海洋环境探测工程技术研究中心，南京 210044；3.国家海洋信息中心，天津 300171)

0 引言

浒苔(Ulvaprolifera)作为绿藻门石莼目的一种大型藻类植物，主要生存在半封闭海域的岩石、潮带滩涂或石砾上。浒苔本身无毒，但当其异常增殖后会影响水中生物生存，恶化水体环境，对海洋渔业及滨海旅游业造成强烈冲击[1-2]。近几年来，浒苔绿潮灾害频发，且主要分布在南黄海的近岸海域。据2018年中国海洋灾害公报[3]统计，2016年、2017年、2018年的黄海沿岸浒苔绿潮最大分布面积依次为57 500 km2、29 522 km2、38 046 km2。2008年的浒苔灾害给当地帆船训练造成了困难，直接影响了青岛奥帆赛的顺利举行，此外还给当地造成了巨大的经济损失，海滩清理、藻类运输和水产养殖损失的经济成本超过2 亿美元[4]。2008年之后，黄海区域浒苔绿潮灾害呈现连续性和常态化的爆发态势[5]。因此，有效地监测预报浒苔已成为控制和减少灾害损失的必由之路。

由于浒苔灾害爆发具有发展迅速、分布面积大等特点，具备宏观性、经济性、多时相性和大范围同步等优势特点的卫星遥感技术成为目前监测浒苔的主要手段。其中，光学遥感观测因其覆盖范围广、光谱信息丰富、空间和时间分辨率高，成为对浒苔灾害进行宏观精准监测的主要数据源[6]。国内外学者针对浒苔的遥感监测提出过许多算法，其中大多数算法采用光谱指数法，比较常用的有归一化植被指数法(normalized difference vegetation index，NDVI)、增强型植被指数法(enhanced vegetation index，EVI)、漂浮藻类指数法(floating algae index，FAI)、替代型漂浮藻类指数法(alternative floating algae index，AFAI)等[7-9]。如Hu[10]通过组合红光-近红外-短波红外波段开发了可消除大气影响的FAI指数，从而利用恰当的阈值准确提取了漂浮藻类信息。Son等[11]发现EVI算法也能够很好地进行浒苔遥感监测。丁一等[12]建立了NDVI与浒苔像元丰度关系的浒苔信息提取模型，利用该模型提取了浒苔的覆盖率。2014年，王宁等[13]对比分析了5种常用的光谱指数对浒苔的监测能力，结果表明，NDVI算法对浒苔灾害各个阶段的监测能力均最强。NDVI指数法是目前自然资源部浒苔灾害卫星遥感应急监测系统的业务化算法。

需要指出的是，虽然光谱指数法是目前浒苔提取的主要方法，但受海洋环境以及卫星成像条件等因素的影响，导致浒苔准确提取的阈值往往存在一定差异。如大部分研究者在利用NDVI指数提取浒苔时所采用的阈值各不相同，有些相差甚远，这可能与不同大气状况、浒苔密集程度及悬浮状态等因素有关[14-23]。这使得利用指数法对浒苔进行自动化监测时，其监测精度往往不稳定，甚至出现严重的误差。此外，也有业务化部门采用人工解译的方法确定准确的阈值，进而提高浒苔监测精度，但这严重影响了浒苔监测的效率，与采用遥感监测的初衷严重不符。因此，开发出可靠的自适应阈值算法，才能从根本上实现指数法对浒苔的全自动监测。

针对以上问题，本文以NDVI光谱指数为例，以GOCI(geostationary ocean color imager)卫星传感器数据为主，提出一种新的光谱指数自适应阈值提取算法，并在浒苔爆发的不同年份以及浒苔生消过程，对提出的自适应阈值算法的适用性进行检验，同时也使用了MODIS(moderate resolution imaging spectroradiometer)数据评价了该算法的适用性，为浒苔的全自动监测提供可靠的技术支持。

1 研究区域与数据算法

1.1 研究区域简介

浒苔主要分布在南黄海海域，因此，以南黄海海域(山东省-江苏省沿海海域)为研究区，南起苏沪分界线，北至山东半岛南岸，经纬度范围为32.5°N～37.1°N，119°E～123°E。胶东半岛成山角到朝鲜的长山串之间的连线将黄海分为北黄海和南黄海2个部分，其中，南黄海是一个近似南北向的半封闭海域，面积约为40.9×104km2。

1.2 卫星数据源和数据预处理

本文从GOCI数据官网下载了L1B数据(http：//kosc.kiost.ac.kr/eng/p10/kosc_p11.html)，利用GDPS这一款针对GOCI数据处理的专业软件[24]，使用了Generate L2/L3 Data工具，将下载好的L1B格式数据处理成了L2数据(主要包括离水辐亮度、遥感反射率、瑞利校正后的遥感反射率等产品[25])，使用瑞利校正后的遥感反射率用于相关研究。

本文也使用了MODIS卫星遥感数据对所开发的自适应阈值算法进行检验。从MODIS数据官网(https：//ladsweb.modaps.eosdis.nasa.gov/search/)下载了L1B数据，对其进行几何校正和大气校正后用于本文的分析。

1.3 监测精度评价算法

在遥感图像提取结果精度评价中，常用的精度评价指标主要有混淆矩阵、总体分类精度(OA)、卡帕系数(Kappa)系数、错分误差、漏分误差等。其中，OA和Kappa应用最为广泛。总体分类精度指被正确分类的类别像元数与总的类别个数的比值，用于评价监测的精度。本研究以目视解译为真实值，使用总体分类精度来评价自适应阈值法对浒苔的提取精度。

2 自适应阈值提取算法开发

本文以常用的NDVI指数为例进行自适应阈值提取算法的开发，其具体表达如式(1)所示。

(1)

式中：RNIR为近红外波段的反射率；RRED为红光波段的反射率。

选取2016年6月25日浒苔爆发时期的GOCI影像，进行自适应阈值算法的开发。首先，对GOCI遥感图像进行了云掩膜、云影掩膜等预处理[26]，进而计算出整幅图的NDVI值，统计得到了如图1(a)所示的NDVI直方图。直方图横坐标为NDVI数值，纵坐标表示整幅图中对应该数值的像元个数。将直方图采用一次、二次等不同次数函数进行拟合，经对比分析后发现，多项式六次拟合函数得到的曲线最接近于原始NDVI直方图，因此，本文选用多项式六次拟合函数，得到的拟合曲线如图1(b)所示，表示如式(2)所示。

f(x)=a1x6+a2x5+a3x4+a4x3+a5x2+a6x

(2)

式中：输出参数a1、a2、a3、a4、a5、a6为六次拟合多项式系数，随着不同的浒苔分布情况而改变。通过拟合曲线图可以看出：浒苔和背景水体有着明显的区分，整个图像在左侧存在一个明显的高峰，经目视解译判别该高峰为背景水体；在右侧存在一块小幅度波动区域，经目视解译判别该波动区域为浒苔。因此，如何确定划分左右两边的分界线，从而确定用于提取浒苔的阈值则是需要解决的关键问题。

为此，本文采用逐步逼近寻找法。首先，通过寻找曲线上的极值，可由其一阶和二阶导数来确定，若记拟合曲线图1(b)函数为f(x)，极值表示如式(3)、式(4)所示。

f′(x0)=0

(3)

f″(x0)≠0

(4)

式中：f′(x0)为一阶导函数；f″(x0)为二阶导函数。在极值点中，若f″(x0)<0，则曲线函数f(x)在x0取得极大值；若f″(x0)>0，则曲线函数f(x)在x0取得极小值。利用该方法找到所有极值点后，从左往右选取曲线的第一个极小值点，该点所对应的横坐标(图1虚线所对应处)即为所求的浒苔提取的NDVI阈值。

图1 自适应阈值法确定的2016年6月25日的浒苔提取阈值

3 自适应阈值提取算法的验证与应用

3.1 自适应阈值提取算法的验证

采用本文建立的自适应阈值提取法，对2016年6月25日GOCI的NDVI指数进行阈值确定，确定的阈值为0.12，再利用该阈值进行浒苔分布提取，结果如图2(c)所示。可以看出，浒苔在爆发时期数量多、分布范围广，这也进一步说明了实现浒苔全自动遥感监测的必要性和迫切性。为了检验自适应阈值提取算法的准确性，对GOCI卫星影像进行了假彩色图像合成，并经目视解译，其中红色像元判断为浒苔分布像元(图2(a))；与此同时，利用固定阈值法也提取了浒苔分布(图2(b))。

由于大部分研究者使用NDVI法提取浒苔信息时采用的阈值各不相同，故本文结合已有文献所使用的NDVI阈值范围以及目视解译结果，将提取浒苔的NDVI固定阈值设为0.15。对比分析以上2种方法和本文自适应阈值法提取的浒苔分布可以看出，自适应阈值法能够很好地将浒苔分布像元进行提取，其提取范围和目视解译以及固定阈值法提取的浒苔分布基本一致，但细致比较后发现与固定阈值法提取结果存在一定的差别，固定阈值方法提取的浒苔面积为4 914.8 km2，而自适应阈值方法提取的浒苔面积则为5 813.5 km2，利用自适应阈值算法提取的浒苔分布更大一些。

为了进一步验证自适应阈值法及固定阈值法提取的浒苔分布，本文以目视解译为真实值，对提取的浒苔分布利用混淆矩阵中的总分类精度进行了精度评价，结果发现固定阈值法提取的总体精度为90.6%，而利用本文所提出的自适应阈值法提取的总体精度为94.7%(表1)。结合以上2个方面的验证得出，自适应阈值法较以往的固定阈值法对浒苔分布的提取精确度有了一定的提高。

图2 2016年6月25日GOCI数据对自适应阈值法的验证

表1 利用固定阈值法和本文自适应阈值法提取的浒苔分布总体精度 %

3.2 自适应阈值提取算法的适用性和稳定性

接下来，为确定自适应阈值提取算法的适用性和稳定性，除2016年6月25日GOCI影像外，进一步选取了2幅不同年份浒苔爆发时期的GOCI影像，时间分别为2017年6月27日和2018年6月24日。运用本文提出的自适应阈值法对浒苔分布进行提取。首先，分析了每一幅影像的NDVI指数直方分布图，并利用本文提出的自适应阈值法确定了各个最优阈值，结果如图3所示。可以看出，2017年6月27日、2018年6月24日2期影像对应的浒苔提取阈值分别为0.22和0.18，与2016年6月25日影像对应的阈值0.12有较大的差别。

图3 不同年份浒苔爆发时期对自适应阈值法的验证

得到不同年份浒苔爆发时期的浒苔提取阈值后，提取出了南黄海海域的浒苔分布图，且将所提取的不同年份浒苔分布进行精度检验。首先，对2017年6月27日和2018年6月24日的GOCI卫星影像进行了假彩色图像合成(图4(a)、图4(d))，同时利用固定阈值法提取了浒苔分布(图4(b)、图4(e))，对比以上2种方法和本文自适应阈值法提取的浒苔分布(图4(c)、图4(f))以及计算出各提取结果的总体精度(表1)，并结合针对2016年6月25日GOCI影像的研究结果(图1、图2)可以发现，在不同年份浒苔爆发时期，运用本文自适应阈值法确定的浒苔提取的NDVI阈值有明显的变化，这些变化势必会造成传统固定阈值方法对浒苔提取结果的误差。由表1可知，利用固定阈值法提取的2016年6月25日、2017年6月27日、2018年6月24日的浒苔总体精度分别为90.6%、86.7%、92.7%，而利用本文自适应阈值法的提取精度分别提升至94.7%、93.4%、95.9%。由此可以得出结论，本文提出的自适应阈值算法对于不同年份浒苔爆发时期的浒苔遥感提取具有良好的适用性。

另外，从提取的浒苔分布结果来看，2016年浒苔爆发情况与2017、2018年相比更为严重，其爆发区域基本涵盖了整个江苏-山东近海区域，而2017年浒苔的情况与之恰恰相反，浒苔爆发现象并不是十分明显，其覆盖的区域较小，主要集中在青岛近海附近。

上文讨论了自适应阈值法在不同年份浒苔爆发期的应用效果，发现其比传统固定阈值法的提取精度更高。接下来将探究自适应阈值法在浒苔生消过程监测中的适用性。选取浒苔灾害的3个典型阶段：初期、爆发、消亡。基于2016年的GOCI数据，选取了浒苔分布由少到多的3期影像，时间分别为5月25日、6月17日、7月16日。运用自适应阈值法，得到各个阶段浒苔提取的NDVI阈值(图5)，2016年5月25日、2016年6月17日、2016年7月16日3期影像对应的浒苔提取的NDVI阈值分别为0.14、0.12、0.27。

图4 不同年份浒苔爆发时期对自适应阈值法的验证

得到浒苔生消各个过程提取的阈值后，提取了南黄海的浒苔分布。为了检验自适应阈值提取算法的准确性，对GOCI卫星影像进行了假彩色图像合成，同时，利用固定阈值法也提取了浒苔分布，对比分析以上2种方法和本文自适应阈值法提取的浒苔分布，结果如图7和表1所示。可以看出，自适应阈值法在浒苔生消过程同样也具有很好的提取精度，尤其是2016年6月17日的数据，自适应阈值提取算法得到的阈值与固定阈值相差较大，对其浒苔提取精度分析发现，固定阈值法的浒苔提取总体精度只有74.4%，而自适应阈值法的总体精度为96.2%，二者相差21.8%之多。以上结果表明，自适应阈值法也适用于浒苔生消过程的监测。

从3期影像提取的浒苔分布图可以看出：在2016年5月25日浒苔发展初期阶段(图6(a)、图6(b)、图6(c))，浒苔在苏北浅滩有零星分布；而在2016年6月17日浒苔爆发阶段(图6(d)、图6(e)、图6(f))，浒苔大规模增殖，在整个江苏近海区域爆发；在2016年7月16日浒苔消亡阶段(图6(g)、图6(h)、图6(i))，浒苔主要在青岛近海海域有较小范围的分布。

由上文的分析可知，对于GOCI数据，本文提出的自适应阈值提取算法具有很好的适用性和稳定性，接下来将开发出的算法应用于MODIS数据，对算法的表现进行检验评价。本文选取了2019年浒苔爆发时期的2景影像，成像时间分别为2019年的6月11日和6月23日。运用本文提出的自适应阈值算法，得到用于浒苔提取的NDVI阈值分别为0.16和0.12。基于此，提取出了南黄海海域的浒苔分布，并对提取精度进行了检验评价，结果如图7所示所示。可以看出，自适应阈值法得到的浒苔分布与假彩色合成图显示的浒苔分布基本一致，利用目视解译对浒苔提取精度的定量分析发现，自适应阈值法的总体精度分别为96.0%和95.4%，这些结果表明本文的自适应阈值法同样也适用于MODIS数据。

图5 自适应阈值法在浒苔生消过程监测中的适用性评价

综上所述，本文提出的自适应阈值法不仅适用于不同年份浒苔爆发时期对浒苔分布的遥感提取，而且也适用于对浒苔生消过程的监测，且对MODIS卫星遥感数据也表现出良好的提取精度。这些结果表明本文的自适应阈值提取算法具有良好的稳定性和适应性。

针对浒苔的遥感提取，此前也有学者提出过自适应阈值的方法。如王蕊等[27]基于数据挖掘提出了浒苔的自适应阈值分区自动检测方法，该方法首先对图像进行分割，再对云覆盖类型进行识别，然后针对不同云覆盖类型，分析反射率与浒苔提取阈值之间的经验关系，进而建立了浒苔的自适应阈值分区检测方法。此外，白羽[28]对比分析了图像分类方法和图像阈值分割方法，进而将图像分类方法和迭代阈值分割法相结合，提出了自适应阈值分割方法，实现了对青岛海域浒苔分布信息的自动提取。可以看出，本文提出的自适应阈值算法与其他学者算法的初衷一致，也即实现对浒苔分布信息的全自动提取，但这些算法的开发思路和所采用的方法均不同；相比之下，本文的自适应阈值算法操作可能更为简单，然而针对不同方法的详细对比分析，仍需要在将来的工作中做进一步的深入研究。

图6 自适应阈值法在浒苔生消过程监测中的适用性评价

图7 自适应阈值法基于MODIS数据适用性评价

4 结束语

本文针对南黄海海域的浒苔灾害，提出了一种新的自适应阈值提取算法，该算法以归一化植被指数NDVI为例，从NDVI直方图分布特征入手，采用曲线拟合、自动查找极值的方法来确定浒苔提取的最优阈值。该自适应阈值提取算法对不同年份浒苔爆发期的浒苔分布提取，以及浒苔生消过程(初期、爆发、消亡)的监测均具有很好的适用性，所提取的浒苔分布较传统固定阈值法提取精度均有一定的提高。

本文所提出的自适应阈值算法可自动确定浒苔提取的最优阈值，提高了浒苔监测的精度，减少了人工目视判别固定阈值的工作量，为实现浒苔的全自动化遥感监测提供了可靠的技术支持。此外，虽然本文提出的自适应阈值提取算法是以NDVI指数为例，仅用于了GOCI和MODIS数据，但从原理上来讲，该方法也能够用于其他光谱指数(如FAI、EVI等)和其他卫星传感器数据(如Landsat、高分等)，这些工作有待进一步开展。