孟苗苗, 郑向阳, 邢前国, 刘海龙
1. 中国科学院海岸带环境过程与生态修复重点实验室(烟台海岸带研究所), 山东 烟台 264003;
2. 中国科学院海洋大科学研究中心, 山东 青岛 266071;
3. 中国科学院大学, 北京 100049
近二十年来, 由大型藻类暴发而导致的海洋灾害在世界范围内呈增长趋势(Morand et al, 2005; Ye et al, 2011; Smetacek et al, 2013), 不仅危害海洋生态环境(Valiela et al, 1997), 还影响海上交通、沿海地区的观光旅游等活动(邢前国 等, 2011)。在中国,大型藻类灾害受到重视始自2008 年夏季, 时值奥运会帆船比赛前, 青岛沿岸暴发了以漂浮大型绿藻浒苔(U. prolifera)为主的绿潮灾害。国外, 类似的海洋大型藻类灾害也在增加, 如法国 Brittany 的绿藻(Schreyers et al, 2021)、墨西哥湾与大西洋的漂浮马尾藻带(Wang et al, 2019; Gower et al, 2019)。
在我国近海与海岸带地区, 黄海夏季由紫菜养殖导致的大规模漂浮浒苔(Liu et al, 2009; Xing et al,2019)以及冬春季黄海、东海的漂浮马尾藻铜藻(Xing et al, 2017; Qi et al, 2017)因规模较大受到了较多关注; 以底栖为主的大型藻也常在本地局部地区形成灾害或环境污染(Xing et al, 2016a; Han et al,2019), 由于在沿海分布广泛, 其生态环境效应与资源管理问题不容忽视。
孔石莼(U. pertusa)是一种潮间带与浅海常见的绿色大型藻, 其出现与水体的高度富营养化密切相关。孔石莼在浅海以底栖为主, 能在非附着的条件下生长, 附着气泡较多时亦可在水面呈漂浮状。孔石莼生长迅速, 生物量大, 常因暴发、过度聚集而形成生态灾害或污染。孔石莼也有着正面的生态功能与资源价值, 在生长过程中可吸收大量营养盐净化水体、吸收二氧化碳, 也可作为食材、动物饲料, 或者深加工作为药物、保健品等原材料(刘岩 等,2010)。
使用卫星影像能够大范围、同步监测绿藻在海面上的分布状态, 为绿潮的防控提供先决条件。微波影像可以根据粗糙度信息识别出漂浮在海面的绿藻(蒋兴伟 等, 2009; Shen et al, 2014), 而热红外影像可以根据不同的热发射率将绿藻从海水中识别出来(Xing et al, 2015)。监测绿藻最常用的是光学影像, 因为绿藻与陆上植被具有相似的光谱特性, 可见光波段反射率低, 红光波段和蓝光波段出现吸收谷, 近红外波段出现反射峰(樊彦国 等, 2015), 具有典型的“红边效应”(安德玉 等, 2018)。早期识别绿藻的遥感算法有差值植被指数(difference vegetation index, DVI)、 比值植被指数(ratio vegetation index, RVI) 、 归一化植被指数(normalized difference vegetation index,NDVI)(Cui et al, 2012; Xing et al, 2018), 其中NDVI 在光学影像绿藻提取中应用最为广泛(Li et al, 2018b)。由于植被指数具有一定局限性, 又有学者提出适用于绿藻提取的指数, 例如漂浮藻类指数(floating algae index, FAI)(Hu, 2009; Wang et al,2016)和虚拟基线漂浮藻类指数(virtual-baseline floating macroalgae height,VB-FAH)(Xing et al,2016a), 这些指数对太阳光、气溶胶变化等的影响敏感性较低, 更适用于海面绿藻的提取。
由于海面绿藻斑块的大小可以从几厘米到几公里不等(Li et al, 2018a), 传统的阈值提取方法属于“硬分类”, 由于混合像元的影响, 在估算绿藻面积时存在一定误差。针对混合像元分解的绿潮面积精细化提取, 有学者基于同步的高分辨率和低空间分辨率卫星影像间分析了混合像元分解(丁一 等,2015), 还有学者通过实测绿藻光谱, 分析绿藻提取指数对绿藻表面覆盖度的响应, 并建立混合像元绿藻覆盖度估算模型(Qi et al, 2016; Li et al, 2018b)。
本文针对孔石莼展开研究, 利用实测同步厘米级分辨率无人机数据提取的绿藻面积, 建立基于Landsat-8 数据计算的植被指数(NGRDI、NDVI和 VB-FAH)以及单波段反射率对绿藻亚像素覆盖度的多个反演模型, 并对模型进行了验证与评价。
1.1.1 卫星与无人机数据采集与处理
本文研究区域位于渤海湾沿岸, 渤海为陆表海,平均水深18m (Xing et al, 2016b), 研究区域现场照片如图1a, 整个研究区海面分布有大量孔石莼, 孔石莼在海里浅水区呈现沉水非定生状态, 退潮时呈干出状态(图1b、c)。本文采用准同步的无人机影像和Landsat TM 卫星影像开展研究, 无人机影像采集时间为2020 年10 月23 日11 时8 分, Landsat 卫星过境时间为10 时49 分。
无人机正射图像采集使用大疆御Mavic 2 为搭载平台, 其搭载的传感器为1/2.3 英寸CMOS, 有R、G、B 三个波段。使用Pix4Dmapper 软件对自带GPS定位信息的无人机影像进行快速校正、拼接处理,对于无法拼接的无特征点图像进行手动拼接处理,输出满足遥感解译需求的影像产品, 拼接后的无人机影像如图3c。对Landsat-8 OLI 数据首先进行辐射定标处理, 再使用FLAASH 大气校正模型以获取反射率数据, 最后进行裁剪处理。为保证地理配准的精度, 将处理后的无人机影像与Landsat-8 影像进行比较, 选取 15 个控制点, 发现配准结果较好,RMSE 为0.51。
1.1.2 高光谱数据采集与处理
为了解孔石莼的光谱特征, 帮助选择参与建模的指数, 于2020 年11 月11 日测量了研究区域不同状态下孔石莼的反射率。采用的光谱仪为 ASD HandHeld2, 波长范围为325~1075nm, 光谱分辨率为1nm; 先测参考板的辐射亮度值, 再测海水与不同状态下孔石莼的辐射亮度值, 每组重复测量3 次取平均值, 反射率计算公式为
式中R( )λ代表被测物体的反射率,L( )λ代表地物辐射亮度值,LP( )λ代表参考板辐射亮度值,ρ( )λ是参考板的反射率。
所测得的海水与孔石莼的反射光谱曲线如图2所示。绿藻的光谱特征为在可见光波段反射率低,在近红外波段反射率高, 与海水光谱的差异明显。将绿藻光谱按照Landsat-8 OLI 传感器第二到四波段的光谱响应函数进行重采样处理, 发现其在蓝、绿、红、近红外波段的反射率低于海水, 在近红外波段反射率高于海水。
1.2.1 无人机影像绿藻提取
在无人机影像上选取一块大小为417×435 像素的样本区, 使用动态阈值法对样本区进行绿藻提取,以人工目视解译结果为基准, 对比提取结果, 进行精度验证。
由于无人机影像只有红、绿、蓝三个波段, 因此选取4 种典型植被指数: 归一化绿红差异指数(normalized green-red difference index, NGRDI)、归一化绿蓝差异指数(normalized green-blue difference index,NGBDI)、红绿比值指数(red green ratio index,RGRI)、过绿指数(excess green, ExG)做绿藻提取研究, 其计算公式如下:
式中R、G、B分别表示红波段、绿波段和蓝波段像素值。
1.2.2 卫星影像绿藻覆盖面积与植被指数关系建模
(1) 植被指数构建
将无人机影像上绿藻的提取结果以30m×30m的矩形为单位建立网格, 使之对应 Landsat-8 卫星影像的每个像素。统计每个像素内绿藻覆盖面积除以该像素的总面积(900m2), 作为每个混合像元的大型绿藻覆盖度(portion of macroalgae, POM)。根据研究区特点选取POM接近1 的像素作为典型绿藻像素,POM 接近0 的像素作为典型海水像素, 端元所在位置如图3a 彩色标注的像素(绿色表示典型的绿藻像素, 红色表示典型的海水像素)。分析其光谱曲线(图4)发现, 在可见光波段, 海水的反射率值高于绿藻, 在近红外波段, 海水反射率值降低, 而绿藻反射率值增高。有明显区分特征的波段是蓝、绿、红、近红外波段。
建立Landsat 卫星影像混合像元POM 与光谱值之间的线性关系:y=ax+b,y为POM,x选取7 种指数参与建模分析, 分别为蓝波段、绿波段、红波段、近红外波段的反射率与NGRDI、NDVI、VB-FAH,其计算公式如下:
式中R为反射率, 下标NIR、RED、GREEN 分别表示近红外、红、绿波段。
(2) 模型建立与验证
在卫星图像上选取两块同样大小的15×12 个像素组成的区域作为样本区与验证区, 如图3a, 为避免配准时的精度偏差, 以3×3 个像素为窗口, 计算每个窗口内的平均POM 与每个窗口内包含的9 个像素的平均光谱反射率,在样本区域建立 POM 与NGRDI、NDVI、VB-FAH、蓝波段反射率、绿波段反射率、红波段反射率、近红波段反射率的回归模型, 选取有较好线性关系的模型公式在验证区域进行验证, 对反演模型拟合的POM 结果与无人机影像提取结果拟合的POM 进行对比分析, 选取最好的反演模型应用到整片无人机覆盖区域。
在无人机影像上均匀选取一定数量的绿藻与海水样本点, 其地表像元亮度值(digital number, DN)均值曲线如图5a 所示。发现绿藻在红波段的DN 值明显低于海水, 因此使用动态阈值法提取绿藻时选用红波段DN 值与4 种植被指数。为了确定基于此研究区域的无人机影像提取绿藻的最佳方法, 选取一块样本区, 以目视解译的结果作为绿藻提取的真实结果, 对比分析5 种指数的动态阈值法提取的结果。
不同提取方法的提取结果如图5b 所示, 从图中可以看出, 使用ExG 指数提取的结果噪点最多。从表1 可以看出, 监测精度最高的是红波段DN 值, 提取一致性达到94%, Kappa 为0.88, F1 score 为0.92。且红波段灰度图(图3d)可以反映出无人机影像单个像素的POM, 为便于之后进行更精细的提取, 使用红波段DN 值对无人机影像进行绿藻提取, 最终提取结果如图3e 所示。
表1 不同绿藻提取方法精度对比Tab. 1 Comparation on accuracy of green macroalgae extraction with different methods
样本区域的建模结果表明, POM 与可见光三个波段的反射率线性关系最好, 其R2都大于0.9, 与近红外波段的线性关系最差。3 个指数中, POM 与NDVI、VB-FAH 的线性关系较好, 与NGRDI 的关系较差(图6)。
选用可见光三个波段所建立的模型在验证区域进行验证。以模型的反演结果作为估测的POM, 以无人机影像提取的结果作为真实的POM, 分析拟合结果(图7)。从表2 可知, 使用绿波段进行反演,R2最大, 为0.92; 均方根误差最小, 平均相对误差也最小, 表明使用绿波段所建立的模型效果最好。
表2 蓝、绿、红波段反演模型结果比较Tab. 2 Comparison of retrieval results for blue, green and red bands’ retrieval models
对Landsat 卫星影像的陆地部分进行掩膜处理, 将绿波段反演公式应用于研究区域, 计算每个像素的POM, 卫星影像反演的绿藻覆盖度如图8a 所示, 无人机影像提取结果反演的绿藻覆盖度如图8b 所示。
计算得到研究区域的卫星反演绿藻覆盖面积为0.4883km2, 无人机影像提取的真实绿藻覆盖面积为0.5729km2, 反演面积与真实面积的比值为0.85。反演结果虽然对部分像素绿藻覆盖度有所低估, 但对于整体绿藻的分布有较好的反映, 可以应用于此区域绿藻覆盖面积的估测。
结合无人机低空航拍影像和Landsat-8 卫星影像, 提出了基于NGRDI、NDVI、VB-FAH 和蓝、绿、红、近红外波段反射率的绿藻亚像素覆盖度估算线性模型。
在参与新模型建模的众多参数中, 蓝、绿、红波段反射率与绿藻亚像素覆盖度呈现出更好的线性关系, 其R2分别为0.95、0.96、0.93。在验证区域对基于此三种模型拟合的POM 结果进行验证,其中绿波段反射率所建立的模型拟合结果最好,其R2为0.92, 均方根误差为0.07, 平均相对误差为10.85%。使用绿波段反射率所建立的模型应用于整个研究区域所反演的绿藻覆盖面积与真实绿藻覆盖面积(无人机图像所提取的绿藻)比值为0.85。此模型可以较好地反映此研究区域的绿藻分布情况, 而对于更广泛的海水区域, 则需要进一步研究。