基于卫星船位数据的北太平洋作业渔船分布及类型研究

孙永文, 张胜茂, 唐峰华, 王书献, 樊伟, 范秀梅, 杨胜龙

（1.大连海洋大学航海与船舶工程学院, 辽宁大连 116023；2.中国水产科学研究院东海水产研究所, 农业农村部东海渔业资源开发利用重点实验室, 上海 200090）

随着卫星技术飞速进步, 获取的船舶实时位置、日期、航速、航向等主要信息愈来愈精准。船位卫星数据被广泛用于船舶安全生产和海洋渔业方面, 如船桥碰撞[1]、船舶碰撞危险区[2]预警及计算渔场捕捞努力量[3]、捕捞强度[4]、渔场分布[5-6]。北太平洋海域辽阔, 但多数渔获物均在特定海域生活, 通过分析该区域的渔船出海作业情况及捕捞上岸的渔获物组成[7-8]判定该区域内主要有围网渔船、灯光鱿钓渔船等[8], 但对于区分不同类型渔船的作业分布的研究较缺乏。为确定北太平洋公海海域不同渔船作业分布及渔船类型, 本文基于2019年卫星数据与实际生产数据, 运用热点分析方法将北太平洋特定渔船作业分布可视化, 再将实际生产数据中的渔船船位与其叠加, 最终通过分析重叠情况初步判别渔船分布与渔船类型。

研究不同类型渔船作业分布情况可为研究不同渔船分布、预测渔场重心提供数据, 有利于掌握研究区域渔场资源分布。目前在渔业应用方面利用热点分析方法研究渔场分布[9-10]、生物多样性[11-14]、渔获物动态变化[15-16]、珍稀海洋生物分布[17]等较普遍。李佳佳等[10]研究太平洋10年涛动（Pacific decadal oscillation index,PDOI）冷期与暖期气候模态下西北太平洋柔鱼渔场的时空分布变化, 对柔鱼资源的可持续开发具有一定科学意义；Barbosa-Filho等[11]研究南大西洋热带海洋生物多样性分布情况, 为保护生物多样性提供数据参考；刘禹希[15]探究了大眼鲷资源分布的时空分布变化,有利于提高大眼鲷资源评估工作的准确性, 为大眼鲷资源的可持续利用及科学管理提供了参考；Bass[17]研究识别海牛死亡率高值的重要空间集群, 为保护海牛提供参考。综上, 利用热点分析方法分析渔船作业分布具有可行性, 结果对管理作业渔船有积极意义。

目前中国海洋一号、海洋二号等多颗卫星搭载数据接收机[18], 北斗三号已实现全面组网并在全球推广使用, 中国自主研发的卫星网络精度提升[19-20], 船位监测能力增强。本文基于卫星获取的渔船船位数据, 经过预处理获得船位记录数量, 同时利用热点分析方法研究区域渔船作业动态变化, 结合渔船实际作业分布情况初步判断渔船类型, 以期为渔业管理提供北太平洋公海海域渔船作业资料, 为未来在该区域内提高捕捞效率以及加强渔业资源的科学管理提供参考。

1 材料与方法

1.1 数据来源

研究区域位于北太平洋公海海域30°—45°N、150°—170°E和40°—45°N、170°—180°E。由于研究区域内渔汛期从每年的3月开始至12月结束[4,7-8,21-22], 本文选择2019年3—12月的数据开展研究。船位数据来源于国家卫星海洋应用中心（https://osdds.nsoas.org.cn/）, 船位记录共14 852 241条, 数据具有27种标准报文, 报文的主要字段包括日期、时间、经度、纬度、水上移动业务标识码（maritime mobile service identifier,MMSI）、航速、航向等。

2019年渔汛期的鱿钓渔船和围网渔船生产数据来源于中国远洋渔业协会渔业专业委员会收集的渔捞日志。鱿钓渔船共有5 905条数据, 围网渔船共有4 109条数据。内容包括作业日期、作业经纬度、日产量等。

根据北太平洋公海渔场的生产资料及研究进展[8,10,23-30]汇总渔船作业范围与作业类型（表1）。目前在研究区域主要有3种作业方式：秋刀鱼渔船采用舷提网网捕方式, 其作业范围基本和鱿钓渔船一致, 是灯光鱿钓渔船的兼捕作业类型, 本研究将舷提网归入鱿钓作业类型进行分析；围网渔船主要作业范围为35°—45°N、150°—170°E；鱿钓渔船采用灯光诱捕方式, 主要作业范围为38°—42°N、152°—179°E。

表1 渔船作业范围与作业类型Table 1 Operating scope and types of light fishing vessels

1.2 数据处理

卫星数据具有27种船位数据的消息格式, 国家卫星海洋中心每天记录的数据中后缀为l1b的文件包含船位数据, 运用python语言完成解析, 解析过程：确定文件类型, 读取需要字段到csv文件中。筛选出本文需要的字段：MMSI、日期、时间、经度、纬度、航速。解析后的数据直接导入SQL Server 2017数据库, 运用语句判断、计算、统计, 数据预处理流程如图1所示, 将预处理后的数据导入ArcGIS软件进行热点分析。

图1 数据预处理流程Fig.1 Data preprocessing procedure

为筛选出准确的渔船信息, 对在该区域作业的渔船按照以下条件进行筛选：围网作业时主船保持低航速向前, 辅船以主船为起点绕一周回到主船；鱿钓渔船作业时保持低航速, 用灯光引诱鱿鱼, 根据作业特点本文界定围网渔船、鱿钓渔船等作业航速相似的渔船作业航速为：0～4节[22,28,33], 作业时间为18∶00—6∶00。以1 h为时间单元, 筛选研究区域内符合条件的作业渔船数量, 计算每小时的平均经纬度作为中心坐标, 以中心坐标计算左上角、右下角的经纬度坐标, 进而计算出球面面积[34]；筛选出符合条件的作业渔船数量与面积的比值作为热点分析的输入要素, 比值的意义为单位面积（km2）的平均船舶数量。为保证验证方法准确, 渔船生产数据的筛选与卫星船位数据的处理方法相同, 限制条件包括相同区域与时间。

式中,S为经纬度球面面积；R为地球半径, 本文取6 371 km；λ1、λ2为经度,φ1、φ2为纬度, 单位rad。

1.3 研究方法

1.3.1 热点分析 为了验证船位卫星数据与实际生产数据, 本文采用局部空间自相关分析方法中的热点分析方法计算得出高值要素、低值要素发生聚类的空间位置, 即热点、冷点, 进而反映要素聚集程度（式2）[4]。

式中,G*i表示要素的空间位置；n表示要素总数量；ωi,j是要素j和要素i之间的空间权重, 要素j和要素i不相邻,ωi,j=0, 反之,ωi,j=1；xj是第j个要素的属性值；Xˉ是为所有要素属性的平均值；S为所有要素属性值的标准差。

使用ArcGIS软件中热点分析工具获得Z得分（GiZScore）、P值（GiPValue）[31]。其中P值表示假设空间要素在整个空间范围内为随机分布的概率, 为2个相邻要素的局部总和与全部要素总和的比值；Z得分表示标准差的倍数, 其绝对值越大说明要素分布越密集, 此时P值也越小, 即出现显著特征, 显著特征为热点或冷点。本文选择置信区间±1.96作为Z得分的阈值, 当Z得分小于-1.96时, 说明该要素有超过95%的可能性是冷点；当Z得分大于1.96时, 说明该要素有超过95%的可能性是热点；当Z得分在-1.96～1.96时, 表示该要素为随机分布的概率较大。根据阈值将研究区域划分为热点区域、冷点区域、随机分布区域, 其中热点区域表示某要素与相邻要素均为高值, 表明渔船每平方千米作业记录高值在地理空间上出现聚类现象的位置, 冷点区域则相反, 随机分布区域属性值特征不明显, 本文主要研究热点和冷点分布变化。

1.3.2 重叠率 为确定热点、冷点与围网渔船、鱿钓渔船的对应关系, 渔船实际作业点与卫星数据记录作业船位的重叠率可作为其判断依据。本文将重叠记录数量与实际作业记录数量的比值作为重叠率（式3）。渔船实际作业数据经过时间和区域条件筛选后, 获得预处理数据, 将预处理后的实际船位与卫星数据记录匹配, 当实际船位处于卫星数据构成的0.1°×0.1°区域[35-36]内, 即属于重叠, 反之不属于重叠。运用ArcGIS软件筛选工具, 匹配的过程如图2所示。

图2 计算重叠数Fig.2 Calculating the number of overlaps

2 结果与分析

2.1 卫星记录的航速统计分析

为确定卫星记录的渔船航速数据情况, 对2019年3—12月渔船航速分布进行统计（图3）, 渔船航速属于0～1节记录有25 676条, 2～3节记录有3 302条, 大于5节的记录均处于5 000～13 000条之间。结合图3, 发现渔船航速0～4节, 基本符合渔船航速的正常作业范围, 并将其定义为渔船生产航速。

图3 2019年3—12月不同渔船航速对应记录数量Fig.3 Number of records corresponding to different fishing boat speed from March to December 2019

2.2 作业渔船热点分析

以Z得分作为热点分析指标, 2019年3—12月平均渔船作业船位数量热点分析如图4所示。其中红色热点分布在37°—44°N、154°—167°E, 在40°—42°N、172°—175°E处有较少数热点；蓝色冷点分布较分散, 多数冷点分布在40°—43°N、150°—153°E, 在30°—36°N、150°—170°E处有少数冷点；白色点为随机分布, 主要分布在36°—43°N、154°—174°E。该结果展示了预处理后卫星数据中渔船作业分布的热点和冷点, 其均具有显著的聚集特征, 可单独分析热点和冷点的分布特点。

图4 2019年3—12月北太平洋作业船位数量热点分析Fig.4 Hot spot analysis on the number of vessels operating in the North Pacific from March to December 2019

2.3 作业渔船热点分布

热点分布的变动代表大多数渔船位置变化情况, 图5所示为渔船作业位置动态变化。可以看出, 5月热点（黄色圆点）分布在39°—41°N、160°—168°E；6月（棕色圆点）多数热点分布在39°—44°N、169°—177°E, 在32°—37°N、166°—171°E处有少数热点；7月（青色圆点）多数热点分布在42°—45° N、164°—168° E, 在43°—44° N、171°—174°E处有少数热点；8月（天蓝色圆点）多数热点分布在38°—44°N、158°—163°E, 在39°—42°N、166°—180°E处有较少热点；9月（紫色圆点）热点分布在38°—45°N、159°—166°E；10月（绿色圆点）多数热点分布在36°—43°N、159°—168°E, 在39°—42°N、169°—176°E处有较少热点；11月（红色圆点）热点分布在40°—43°N、154°—157°E。分析各月热点分布变化可知：多数渔船先向东北移动, 之后由东北返回西南。在8月有少数热点抵达最东部, 在7—10月多数渔船分布在研究区域中部, 11月移动至最西部。观察图4发现, 北太平洋作业渔船主要分布在37°—44°N、154°—176°E, 热点作业渔船的动态变化与往年北太平洋鱿钓渔船作业分布及轨迹变化基本相似。

图5 2019年5—11月北太平洋作业船位数量热点分析Fig.5 Hot spot analysis of the number of operating vessel positions in the North Pacific from May to November 2019

2.4 作业渔船冷点分布

冷点分布代表相对少数渔船的分布情况, 如图6所示, 3月（深蓝色圆点）、4月（红色圆点）冷点分布在34°—36°N、157°—159°E；5月（黄色圆点）冷点分布在35°—38°N、160°—164°E；6月（棕色圆点）冷点分布在36°—38°N、155°—165°E, 分布较分散, 出现以上情况的原因应是这几月不断有渔船加入到捕捞队伍中导致；7月（青色圆点）冷点分布在42°N、152°E附近；8月（天蓝色圆点）冷点分布在39°—43°N、154°—155°E；9月（紫色圆点）冷点分布在39°—43°N、153°—155°E；10月（绿色圆点）冷点分布在41°—45°N、150°—156°E；11月（红色圆点）冷点分布在41°—43°N、150°—156°E。观察各月冷点分布可知：3—6月少数渔船整体先向北移动, 7—11月向渔场（41°—44°N、150°—156°E）西北部移动。经分析, 船位数据冷点初期在南部分布, 后期转移至41°—44°N、150°—156°E, 与往年北太平洋围网渔船作业分布位置及变化轨迹基本相同。

图6 2019年3—11月北太平洋作业船位数量冷点分析Fig.6 Cold spot analysis of the number of working positions in the North Pacific from March to November 2019

2.5 作业渔船类型验证

为了进一步验证作业渔船的类型, 经过重叠条件预处理筛选可得热点与鱿钓渔船重叠2 763条作业记录, 冷点与围网渔船重叠795条作业记录。根据作业范围筛选出冷热点数量, 其中卫星渔船数据记录热点4 498条、冷点1 536条。运用Python语言对进行计数, 得到重叠的实际作业数据记录数：围网1 017条、鱿钓3 300条。将实际鱿钓渔船、围网渔船作业船位数据分别与热点（图7A）、冷点叠加（图7B）可视化, 并计算重叠率得：鱿钓渔船为83.7%, 围网渔船为78.1%。由此可知, 卫星数据记录作业船位与实际作业船位有较高的重叠率。覆盖率结果表明, 热点应为鱿钓渔船, 冷点应为围网渔船, 结果得到验证。

图7 热点、冷点与渔船实际作业点叠加Fig.7 Overlay of hot spots, cold spots, and actual operating points of light fishing boats

3 讨论

3.1 作业渔船动态分布分析

通过分析发现, 基于2019年3—12月渔船卫星数据得到的渔船作业船位动态热点分布结果基本符合渔汛期渔船在研究区域动态分布变化, 利用船位卫星数据中经纬度确定渔场分布位置变化, 实为渔船作业船位分布, 说明两者可作对比, 结果可通过渔场重心的变动对照作业渔船船位数据冷热点变化开展辨析。在考虑相同限制区域、航速、作业时间的情况下, 本文研究2019年3—12月北太平洋作业渔船热点分布与利用其他数据及方法研究2019年之前的北太平洋渔场重心变化有较高的相似性。如徐博等[21]、戴澍蔚等[29]基于海洋环境数据与渔捞日志, 利用广义加性模型提出2014—2017年灯光围网渔场重心集中分布在39°—43°N、147°—155°E范围内, 渔场重心先在北太平洋公海南部、后向东北移动, 9月以后逐步向西南移动, 与本研究的冷点分布变化类似。魏广恩等[23]、王韫沛等[25]基于海洋环境数据与渔捞日志, 利用聚类分析、神经网络分析2004—2015年正常气候下39°—45°N、150°—165°E海域的鱿钓渔场重心变化规律, 发现在8—9月向东北移动、9—10月向西南移动, 与本研究的热点分布变化类似, 但也有明显区别, 如渔场重心在9月底向西南方向(39°N、148°E)折返[23], 而作业渔船在10月仍有大量热点分布于东北（40°N、172°—175°E）, 作业渔船冷点分布中10月仍有大量冷点处于东北（43°N上下、156°E左右）。因此可知, 上述鱿钓渔船及灯光围网渔船的渔汛期渔场重心变化结果与本文卫星数据的研究结果相似, 且上述研究的数据全面, 针对作业渔船动态分布分析方法相对成熟等, 故可用于验证本文结果。

3.2 卫星数据与渔船类型分析

基于卫星船位数据对渔船作业分析研究, 再结合海洋环境因素、渔捞日志, 利用合理方法可了解渔场重心具体区域及其变化、不同海洋生物的洄游路线等, 可提高渔场重心的准确率、捕捞效率[21,25,37]。本研究通过计算卫星渔船船位数据与渔船实际作业位置重叠率, 验证了渔船作业类型的判别, 即热点应为鱿钓渔船, 冷点应为围网渔船。研究过程属于在卫星船位数据基础上深入挖掘工作, 并对不同渔船作业方式进行区分判别, 弥补了该方面的空白。北太平洋公海存在3种作业方式, 即灯光鱿钓[38]、灯光围网[39]、舷提网作业[40], 其中有一部分渔船是鱿钓和舷提网兼捕作业, 该类兼捕作业渔船在每年9—11月主要开展舷提网作业[41], 其他月份开展鱿钓作业, 所以基本无法区分这2种渔船, 故研究通过冷热点分析仅区别了围网与鱿钓2种作业类型, 若要区别灯光鱿钓与舷提网作业类型不能仅依据卫星船位数据, 更需要根据作业时间与灯光亮度（鱿钓用传统卤素灯而舷提网用红白色LED灯）等来辅助判别。另外通过上述船位重叠率计算得到有少部分作业船位数据并未重叠, 分析可能与选取航速0～4节有关, 范围偏广, 筛选后仍包含了渔船停船未作业和少量渔船找鱼的数据；渔捞日志提供了重复的实际生产作业船位位置, 记录船位太少。若要精准识别渔船类型与具体渔船作业分布大, 可通过夜光遥感技术手段[41]和每艘渔船对应的MMSI识别辅助验证。

本文基于2019年卫星渔船船位数据得到渔船在研究区域的作业动态分布情况, 再利用中国围网渔船、鱿鱼钓渔船生产数据叠加卫星船位数据的方法, 达到围网渔船、鱿钓渔船不同作业类型的初步识别。本研究利用数据库进行数据预处理到利用ArcGIS软件热点分析整个过程中, 其中作图过程花费时间较长, 后期将从使用程序实现自动化出图方面考虑, 提高出图效率。另外,利用热点分析方法将有效卫星船位数据按置信区间±1.96分为3个区间, 其中随机分布区间（白色点）没有明显规律, 有可能是部分不同作业类型的船只混杂在一起, 导致渔船类型难以判断, 未来研究需要借助更好的技术手段来辅助识别。