基于Landsat8数据的近海养殖区自动提取方法研究

武易天， 陈 甫， 马 勇， 刘建波， 李信鹏

(1.中国科学院遥感与数字地球研究所,北京 100101； 2.中国科学院大学,北京 100049)

0 引言

近海水产养殖作为我国渔业经济的重要组成部分，长期以来保持着高速增长势头。根据联合国粮农组织2014年发布的数据，2012年我国水产品年产量占全球渔业市场的61.69%，其中藻类海产品总量占全球市场的53.97%[1]。然而，随着近海水产养殖业的高速发展和养殖面积的不断扩大，大规模养殖浮绳的散乱分布严重影响了海上交通； 加之人工饲养残留在海水中的饲料会造成近海岸水体的富营养化，严重破坏近海岸区域的生态平衡，造成重大损失。因此，为了维护近海岸区域的生态环境以及当地水产养殖业的可持续发展，促进近海资源的合理利用和综合管理，对近海岸水产养殖区域进行快速、准确的监测是很有必要的[2]。

卫星遥感技术以其周期短、重复观测且能够大面积获取地物空间信息等特点，已成为地表动态监测的重要手段。利用卫星遥感影像能够快速、准确地对近海养殖区域进行动态监测，在近海养殖区域的科学规划与管理中发挥重要作用。马艳娟等[3]基于ASTER卫星数据，通过分析影像内各波段对水产的敏感性差异，构建了归一化水产指数(normalized difference aquaculture index，NDAI)，并结合深海提取指数(marine extraction index，MEI)完成了近海水产养殖区养殖水体与其他水体的区分； 周小成等[4]通过构建卷积算子，使用邻域分析法增强水产养殖区的空间纹理信息，对ASTER影像中的水产养殖区域进行了提取； 林桂兰等[5]利用方差法并结合纹理分析，对厦门海湾的养殖区域进行提取，得到了养殖专题图； 孙晓宇等[6]采用面向对象的方法，综合运用养殖区的光谱、纹理、形状特征及空间邻近关系等信息，并结合多时相遥感数据，对珠江口海岸带地区水产养殖场的变化情况进行了提取。然而，上述方法识别的都是养殖区的外部范围，并未实现对养殖区内部的精确提取； 而养殖区外部往往有许多水域并未进行养殖生产，造成了这些方法的结果都出现了较大误差，很难满足监控管理的实际需求。卢业伟等[7]基于RapidEye高空间分辨率影像，提出了一种基于近海养殖区的特征指数并结合纹理信息的养殖区自动提取方法，对福建三都澳养殖水域进行高精度提取，很好地识别出养殖水域间的非养殖水域，大大提高了提取精度。

但是，由于近海岸水域受悬浮泥沙浓度和叶绿素浓度等的影响，常呈现出不同的反射特性，导致近海养殖区水域背景复杂多变； 加之云雾的影响，单一的水产特征指数在复杂水色背景下很难取得较好的提取效果。为了解决复杂背景下水产养殖区的高精度提取问题，本文首先利用基于正交子空间投影的约束能量最小化(orthogonal subspace projection-weighted constrained energy minimization，OWCEM)算法[8]对水产养殖区域进行增强，消除复杂背景对水产养殖区提取的干扰； 然后提出局部标准差自适应分割法(standard deviation adaptive segmentation，SDAS)对水产养殖区域进行分割； 最后基于分割结果实现近海养殖区的高精度自动提取。

1 研究区概况及数据源

1.1 研究区概况

为了验证算法的适应性，本文选取福建宁德三都澳和山东荣成桑沟湾2处近海水域进行水产养殖区提取实验。这2个区域的近海水产养殖区分布广泛，养殖密度大且水产品种多样，是我国著名的水产养殖模范基地。2个实验区域在气候和水文条件等因素上有一定的差异，适合用于测试水产养殖区提取算法的稳定性和可靠性，其分布位置如图1所示。

(a) 山东桑沟湾养殖区 (b) 福建三都澳养殖区

图1实验区覆盖范围示意图

Fig.1Sketchmapforcoverageofexperimentalareas

山东桑沟湾水产养殖区位于山东半岛东端，除东面湾口与黄海相连外，其余三面被陆地环绕，主要养殖水产为贝壳类(箱式养殖)和海带(浮绳式)。福建三都澳水产养殖区位于浙江省宁德市，是我国典型的近封闭型港湾，湾内水产养殖区域众多，主要养殖种类为紫菜(浮网式)、海带(浮绳式)和鲍鱼贝壳类(网箱式)等。

不同种类的水产养殖方式有一定的差别。如图2(a)所示，鲍鱼和其他贝壳类常放在渔排下悬吊的养殖网箱中，而渔排多由木板和塑料材质搭建，其反射特性与水体差异明显，在遥感影像中最容易提取，因此不是本文探讨的重点。不同于鲍鱼养殖，紫菜在海水中漂浮时间过久便会腐烂，在落潮后需要进行干燥晾晒，所以渔民多利用高度20 m左右的大毛竹结合浮网在沿海滩涂养殖紫菜(图2(b))。而图2(c)和图2(d)反映的则是海带的养殖方式，相比于紫菜，海带在收割前并不需要晾晒，因此多将其系在长绳上进行养殖，而长绳要等间隔地系有漂浮球，使海带不会漂浮到海水较深处，并保证其光合作用的效率。

(a) 鲍鱼养殖区 (b) 紫菜养殖区 (c) 海带养殖区1 (d) 海带养殖区2

图2水产养殖区不同的养殖方式

Fig.2Differentcultivationmodelsofaquacultureareas

1.2 数据源

遥感影像选取美国Landsat8卫星数据，通过美国地理地质调查局全球数据服务网站(http: //earthexplorer.usgs.gov/)下载，空间分辨率为30 m。结合水产的养殖周期，选取水产生长较为茂盛的时间节点进行水产养殖区提取实验，其中山东桑沟湾海域选择2016年4月15日获取的数据，福建三都澳地区选择2014年4月10日获取的数据。而在精度验证中则使用多光谱与全色波段融合后的Landsat8数据(15 m空间分辨率)，以及从中国资源卫星中心获取的2016年06月17日的桑沟湾海域高分二号(GF-2)数据和2014年05月08日的三都澳海域高分一号(GF-1)数据。

2 原理与方法

为实现复杂背景下的水产养殖区提取，本文首先利用OWCEM算法对水产养殖区域进行增强，然后结合纹理信息提出SDAS对水产养殖区域进行分割提取，从而实现近海养殖区的高精度自动提取。OWCEM- SDAS方法的整体流程如图3所示。

图3 水产养殖区提取方法(OWCEM-SDAS)流程图Fig.3 Flow char of extraction method(OWCEM-SDAS) for aquaculture areas

2.1 数据预处理

对获取的Landsat8原始影像进行以下预处理：

1)大气校正。首先对Landsat8影像数据进行辐射定标，然后对其进行FLAASH大气校正，以消除影像数据中辐射亮度的各种失真。

2)几何纠正。几何纠正的目的是通过一定数量的控制点，进行不同影像间的像元插值运算和空间变换，以达到消除遥感影像几何畸变的目的。本文通过在遥感影像中均匀布设控制点，对遥感影像进行几何精纠正。

3)水陆分离与影像裁剪。在以上基础上，根据行政区边界和全球30 m空间分辨率地表覆盖产品[9]对遥感影像进行水陆分离，裁剪出实验区域。

2.2 OWCEM算法养殖区增强

在近海岸水域，海水中叶绿素和悬浮泥沙浓度的分布差异使传统方法所构建的水产特征指数很难得到水产养殖区增强的效果，会显著影响后续的纹理分析和水产养殖区提取过程，而OWCEM算法能够对多维特征空间进行目标信号增强和识别，并广泛应用于高光谱遥感影像目标识别中[10]。

2.2.1 样本选择和波段扩展

对多光谱影像使用OWCEM算法前，需要进行波段扩展[11]。扩展的波段主要来源有2类： ①与目标地物有关的波段特征信息，在水产养殖区提取中主要由对水产较为敏感的波段构成； ②通过计算影像中其他像元与端元样本光谱间的相似性所构成。

为了突出水和水产养殖区间的差别，本文在光谱特征空间内新扩展了3个水体指数，分别是改进的归一化差异水体指数(modification of normalized difference water index，MNDWI)[12]，以及Feyisa等[13]提出的水域自动提取指数(automated water extraction index，AWEI)——AWEInsh和AWEIsh，并增加了叶绿素浓度指数(chlorophyll-a concentration，CHL)和悬浮泥沙浓度指数(total suspended matter concentration，TSM)[14]，用以减少水色差异对水产养殖区提取过程的干扰。在相关性特征空间中，为了计算端元样本和其他像元间的相似性，选取光谱角距离SAD，相关系数corr，光谱距离s和光谱信息散度SID[15]这4个特征指数用于有效区分影像中的像元点与先验样本间的差异，提高水产养殖区提取的精度。上述9个指数的表达式如表1所示，其中，AWEInsh与AWEIsh采用文献[8]的方法进行归一化。

表1 波段扩展中所包含的特征指数Tab.1 Feature indexes included in band expansion

(续表)

表1光谱特征中的band1，…，band7分别为Landsat8多光谱数据中的第1波段，…，第7波段； 相关性特征中的x为先验样本的特征向量，y为影像中除先验样本外的所有像元点特征向量；p和q分别为x和y的概率密度分布； 光谱信息散度中的l为第l波段，l=1,2,…,L,L为影像的波段总数，这里取L=7;xl和yl分别为像元x和y中l波段的值，因为在数据预处理中已进行了辐射和大气校正，因此这里的值代表了地物像元在不同波谱范围上的地表反射率。

完成波段扩展后，得到9个新的特征波段，再加上Landsat8影像原始的7个波段，总计16个波段。这16个波段不仅包含了水产养殖区的特征信息，还包含了非养殖水域、悬浮泥沙和叶绿素浓度的信息，综合分析这些特征信息可以更加容易地区分水产养殖区与非养殖水域，提高水产养殖区提取精度。

2.2.2 OWCEM养殖区增强

经过波段扩展后得到的特征空间并不能够直接进行水产养殖区提取，还要经过OWCEM增强处理后结合自适应阈值进行养殖区分割和提取。OWCEM算法的数学表达式为

(1)

式中：N为像元总数;xi为影像中每一个像元的观测信号，xi=(xi1,xi2,…,xiL)T，1≤i≤N;L为影像波段数(通常情况下L<

(2)

式中P为进行正交子空间投影的系数，即

P=I-dd+,

(3)

式中：I为单位矩阵；d为输入水产样本；d+为d的伪逆。

将R带入式(1)可得

(4)

OWCEM算法的具体求解过程如图4所示。

图4 OWCEM算法流程图Fig.4 Flow chart of OWCEM algorithm

用OWCEM对近海水产养殖区域进行增强的步骤如下：

1)首先，由Landsat8影像选取原始水产样本d，在波段扩展后利用式(2)计算正交投影加权自相关矩阵R。

2)其次，在求解出R后，将其带入到式(3)中即可求得OWCEM能量最小解wOWCEM。

3)最后，将该系数与波段扩展后的影像相乘，便可得到增强后的近海水产养殖区域图。

需说明的是，增强后的结果只是归一化到(0,1)区间上的灰度影像，而水产养殖区域的精确提取还需要进行自适应阈值分割处理。

2.3 SDAS法养殖区提取

在大范围水产养殖区遥感精确提取中，受悬浮泥沙等因素的影响，在全局范围内很多养殖区与非养殖区在光谱特征上仍难以区分，这就需要结合空间信息进行进一步分割提取。由于水产养殖区呈条带状间隔分布特征，因此本文在使用OWCEM算法增强养殖区后，采用SDAS法对养殖区进行自动精确提取。

通过对OWCEM算法增强后的影像进行分析，发现水产养殖区呈固定宽度的条带状分布，其宽度在Landsat8影像中多为3个像元，而且相邻养殖区域间被非养殖水域隔开。基于这种纹理特征，本文采用5像元×5像元的窗口对增强后的影像进行分割，并计算各小窗口内部像元灰度分布的标准差。由于所选窗口宽度大于单条养殖区域，因此各窗口内像元分布只可能有养殖区域与非养殖区交错分布或者该窗口内部不存在养殖区域2种情况。在养殖区与非养殖区交错分布的情况下，分割窗口的标准差较大； 而在不包含养殖区域的情况下，局部窗口标准差则较小。因此，可以使用OWCEM处理中选取的样本所在局部窗口的标准差作为阈值，剔除不包含养殖区的小窗口。对剩余的包含水产养殖区的局部窗口，使用最大类间方差分割法Otsu[16]进行水产养殖区的精确提取。

Otsu法是一种自适应的阈值选择方法，通过使影像目标类和背景类二者的类间方差最大化，自动求得阈值，对影像中的目标和背景进行分割。在影像分割时，假设目标像元占影像总像元比例为w1，且目标像元类别的平均灰度为u1； 对应地，影像中背景像元所占影像总像元个数比例为w2，其类别平均灰度为u2，则影像的总平均灰度u为

u=w1u1+w2u2。

(5)

由此得到目标和背景二者的类别间方差g的数学表达为

g=w1(u1-u)2+w2(u2-u)2,

(6)

进而可求得

g=w1w2(u1-u2)2。

(7)

通过选择分割阈值T，使g最大，即

(8)

式中DNmax为影像的最大灰度值。

在局部标准差滤波后的窗口内自动选取阈值，对养殖区与非养殖区进行分割，可以得到近海养殖区的精确提取结果。相对于仅通过全局阈值提取近海养殖区的面向对象提取法，局部标准差自适应分割法能够充分利用养殖区局部的空间分布特性，根据不同子区域自动求得分割阈值，从而很好地减少了悬浮泥沙和叶绿素浓度对分割阈值选择的干扰，提高了近海水产养殖区的提取精度。

3 实验结果与精度检验

3.1 对比实验的选择

由于现有的提取算法鲜有基于Landsat8数据进行近海水产养殖区精确提取的，为了保证实验的可靠性，本文根据Landsat8的数据特征，参考已有的其他卫星数据的水产养殖区提取方法，在2个实验区进行了面向对象水产养殖区提取实验，并与本文实验进行对比分析。为了进一步提升面向对象方法的提取精度，在设置提取阈值时，本文在预设阈值的左右再进行10步(步长为0.02)搜索，以选取整体提取精度最好的阈值作为提取参数。

3.2 实验结果分析

3.2.1 山东桑沟湾提取结果

图5示出山东桑沟湾水产养殖区提取结果。图5(a)是山东桑沟湾水域Landsat8影像的B4(R)，B3(G)，B2(B)波段合成影像； 图5(b)是图5(a)中所标注的3个局部矩形区域A，B，C的放大显示图； 图5(c)是使用面向对象提取方法对山东桑沟湾海域进行水产养殖区提取后的结果，图中黄色区域为数据预处理中通过掩模运算去掉的陆地区域，黑色区域为非养殖水域，白色条纹状区域是水产养殖区。图5(c)选择了图5(a)中同样的3个矩形框，对局部区域进行放大，如图5(d)所示。而图5(e)则是用本文所提方法在山东桑沟湾水域提取的水产养殖分布结果，图中同样选了3个区域放大(图5(f))，以便进行对比分析。

从图5显示的结果来看，影像中对应的3个区域A，B，C都有一定程度的悬浮泥沙干扰。其中区域A范围较大，受悬浮泥沙影响最大，右侧的部分条状养殖区几乎不可见，仅能看到模糊的轮廓； 区域B中的水产纹理较为密集受叶绿素浓度差异的影响，左下角灰度差异较为明显； 而区域C虽然也有悬浮泥沙的影响，但该区域内悬浮泥沙浓度分布较为均一，对水产养殖区提取的影响相对较少。从目视效果上看，本文方法的提取结果(图5(e)和(f))纹理较清晰，水产养殖区内部间隔的非养殖区水域也能够被有效地剔除； 而面向对象的提取结果(图5(c)和(d))则出现了较高的误判率，其中在区域A中右侧被悬浮泥沙影响的水产养殖区几乎没有提取出来，区域B中的结果中存在大量相连的区域，区域C中因全局阈值的影响，提取出的区域也相对有限。需要说明的是，受Landsat8影像空间分辨率的影响，矩形度等特征指数在面向对象提取过程中所起到的作用相对有限，这也增加了面向对象提取水产养殖区域的难度。

(a) Landsat8影像 (b) (a) 中3个区域放大显示

(c) 面向对象提取结果 (d) (c) 中3个区域放大显示

(e) 本文方法提取结果 (f) (e) 中3个区域放大显示

图5山东桑沟湾水产养殖区提取结果

Fig.5AquacultutreareaextractionresultinSanggouBay,Shandong

3.2.2 福建三都澳提取结果

福建三都澳地区的水产养殖区域提取结果如图6所示。图6中各子图的位置与图5类似，图6(a)，(c)和(e)分别是福建三都澳水域的Landsat8影像、面向对象提取结果和本文方法所提取的结果，图6(b),(d)和(f)分别是图6(a),(c)和(e)中3个矩形区域的放大图。从影像中可以看出，福建三都澳海域受悬浮泥沙浓度的影响要比山东桑沟湾地区的小，但该景影像中有部分薄云的干扰，导致部分区域反射率升高，而左侧的区域B未受云影响，水产和非水产养殖区的整体反射率较低。在提取结果中，面向对象的水产养殖区提取出现了很多错分的情况，其主要原因在于大范围区域内的某些非养殖水域反射率和部分养殖区反射率重叠，因此在设置水产养殖区提取的阈值时很难通过一个全局的阈值对大范围水产养殖区进行提取。而本文方法在使用OWCEM进行目标信息增强后，在分割提取时采用2层不同尺度的提取，从目视效果上看获得了较好的结果。

(a) Landsat8影像 (b) (a) 中3个区域放大显示

(c) 面向对象提取结果 (d) (c) 中3个区域放大显示

(e) 本文方法提取结果 (f) (e) 中3个区域放大显示

图6福建三都澳水产养殖区提取结果

Fig.6AquacultutreareaextractionresultinSanduaoBay,Fujian

比较图5和图6可以看出，与面向对象提取法相比，本文方法较好地利用了近海水产养殖区的光谱-空间信息。但本文方法在图6的区域B中也存在部分错分现象，并非所有养殖水域间的间隔(非养殖水域)都能很好地被提取出来。从整体而言，面向对象方法由于受全局阈值的影响，其整体分类结果较为一般； 而本文提出的分类方法提取的水产养殖区纹理清晰，精度较高，具有较高的实用价值。

3.3 精度验证

验证数据都使用ENVI软件进行预处理。数据与处理算法精度的验证是通过在每个实验区内随机生成一定数量的验证样本，然后结合同一地区的高空间分辨率影像数据(融合后的Landsat8和GF-1/2影像)，对随机生成的验证样本点进行目视解译和分类。在进行目视解译前，对GF数据进行了辐射和大气校正，且与Landsat8数据进行了几何配准。最后利用验证样本点与近海养殖区提取结果构建混淆矩阵，得到提取后的总体分类精度和Kappa系数以及水产养殖区的制图精度和用户精度。验证点的选取和目视解译参考影像如图7所示。图7(a)为山东桑沟湾的近海养殖区域图，图中绿点为非养殖水域的验证样本点，红点为水产养殖区的验证样本点。图7(b)为结合GF-2数据对图7(a)中验证样本点进行目视解译和分类的结果。同理，图7(c)为福建三都澳水域的近海养殖区域图，图中绿点为非养殖水域的验证样本点，红点为水产养殖区的验证样本点，而图7(d)为图7(c)中验证样本点在对应区域GF-2数据上的展示。可以看出该景GF-2数据受水色和云雾影响，虽然可以作为参考数据对Landsat8的验证样本点进行目视解译，但较难直接用作水产养殖区的精确提取。在对2个实验区验证样本点目视解译后，对水产养殖区提取结果计算混淆矩阵，其结果分别如表2和表3所示。

(a) 山东桑沟湾养殖区影像(b) 山东桑沟湾养殖区部分验证样本点

(c) 福建三都澳养殖区影像(d) 福建三都澳养殖区部分验证样本点

图7 验证样本点在实验区参考影像中的分布Fig.7 Vertification samples of test regions distributed in reference images表2 山东桑沟湾区域提取结果混淆矩阵Tab.2 Confusion matrix of aquaculture area extraction result in Sanggou Bay, Shandong

表3 福建三都澳区域提取结果混淆矩阵Tab.3 Confusion matrix of aquaculture area extraction result in Sanduao Bay, Fujian

由表2和表3可以看出，面向对象方法在提取时将大量的非养殖水域像元误提取为水产养殖区，而本文方法误检的像元个数相对较少； 在非养殖水域的识别上，2类方法表现都较好，误检像元较少。

2个实验区2种提取方法的整体提取精度结果如表4所示。

表4 面向对象和本文方法提取精度对比Tab.4 comparison of extraction accuracy using obeject-oritented method and proposed method (%)

从表4可以看出，在山东桑沟湾的水产养殖区提取过程中，面向对象方法用户精度达到88.78%，但其制图精度很低，Kappa系数仅为57.03%，总体精度仅为69.69%； 本文方法提取精度较高，总体精度达到了94.86%。在福建三都澳实验区，面向对象提取精度有所提高，但总体精度仍未达到80%； 而本文方法虽然精度相比山东桑沟湾实验区有所下降，但Kappa系数仍在90%以上，总体精度也比面向对象方法高出15%。

综上所述，在山东桑沟湾和福建三都澳水产养殖区提取实验中，面向对象提取算法总体精度不足80%，而本文算法总体精度维持在93%以上。证明本文方法能够在复杂水色背景下，对水产养殖区进行精确提取。

4 结论

利用遥感技术进行近海养殖区的提取已经广泛地应用于生产实践，但目前传统的提取方法在海域水色较复杂的情况下提取养殖区的精度较低，远远不能满足实际应用的需求。本文在充分利用了近海养殖区的光谱和空间信息的基础上，利用正交子空间投影的约束能量最小化(OWCEM)方法，对复杂海域下的水产养殖区进行增强； 并通过局部标准差自适应分割(SDAS)算法实现了养殖区的自动分割和提取，大大提高了养殖区提取结果的精度。与现有的近海养殖区提取方法相比，本文方法具有以下优点：

1)充分利用了影像中所有的波段信息，能够在复杂水色背景下有效增强近海水产养殖区，降低复杂背景对提取的影响。

2)自动化程度高。除算法初期需要输入少量样本点(不超过10个)外，整个处理流程完全自动进行，大大减少了提取的工作量，也降低了方法使用难度。

3)基于局部标准差自适应分割方法自动调整每个子区域的分割阈值，解决了全局分割阈值不能精确提取水产养殖区的问题。

但本文方法还存在某些不足，在水产养殖区提取过程中未考虑水产生长情况对提取的影响，导致部分区域出现了一些错分的现象，这也是下一步进行深入研究的目标，即结合多时相遥感数据更有效地对近海养殖区进行提取并监控水产的生长情况。