西北太平洋柔鱼丰度的灰色灾变预测

解明阳，陈新军,2,3,4,5*

(1. 上海海洋大学 海洋科学学院，上海 201306；2. 农业农村部大洋渔业开发重点实验室，上海 201306；3. 国家远洋渔业工程技术研究中心，上海 201306；4. 大洋渔业资源可持续开发教育部重点实验室，上海 201306；5. 农业农村部大洋渔业资源环境科学观测实验站，上海 201306)

1 引言

柔鱼（Ommastrephes bartramii）是大洋经济头足类，广泛分布于西北太平洋，是我国鱿钓渔船的重要捕捞对象[1]。自1993年我国开始开发柔鱼资源以来，捕捞规模与作业海域不断扩大，年产量基本稳定在6万～10万t[2]。由于柔鱼是一种“生态机会主义”物种，其中心渔场与资源丰度极易受海洋环境与气候变化的影响[3]，因此柔鱼资源量会出现剧烈的年间变化。已有研究表明，柔鱼的资源量变化与环境因子变化，例如海表温度（Sea Surface Temperature, SST）[4]、黑潮[5]有密切的关系；资源量的骤增骤减现象也受气候变化的影响，例如 El Niño-La Niña 现象[6]、太平洋年代际涛动（Pacific Decadal Oscillation，PDO）现象[7]等影响。目前国内外学者已经对北太平洋柔鱼渔场时空分布[8]、资源评估与管理[9]、资源量预报模型[10]等做了较为全面的研究。对柔鱼资源量进行年度的灾变预测，将有利于掌握西北太平洋柔鱼的资源状况，为资源合理利用提供依据。

灰色系统理论是一门有关不确定系统的理论，其优点在于允许样本数量较少且服从任意分布，这对于缺乏数据的渔业科学领域来说有很好的应用前景[11]。目前，灰色系统理论在资源量预测[12]、渔场空间分布[13]、渔业水质评价[14]等方面已取得了较好的成果。为此，本研究采用灰色系统研究方法中的灰色灾变预测模型，对西北太平洋柔鱼资源丰度的丰歉年进行预测，以达到预警目的，为北太平洋渔业生产企业和管理部门及时进行有效管理提供技术支撑。

2 材料与方法

2.1 渔业生产数据

渔业生产统计数据来自上海海洋大学中国鱿钓技术组，时间为1995-2017年，空间范围为35°～45°N、140°～179°E，统计内容包括日期、经度、纬度、日产量和作业渔船，空间分辨率为1°×1°。以单位捕捞努力量渔获量（Catch Per Unit Effort，CPUE）表征柔鱼资源丰度，以1°×1°为一个渔区，计算各年每个渔区的CPUE值。

SST数据来源于美国国家海洋和大气管理局（National Oceanic and Atmospheric Administration,NOAA）的网站（http://oceanwatch.pifsc.noaa.gov/erddap/index.html），时间范围为 1995-2017 年的1-12月，时间分辨率为月，利用平均法计算各年作业渔场平均SST。

2.2 CPUE标准化

统计每年所得到的总渔获量（Catch），以每年捕捞的船数作为捕捞努力量(Effort)，计算名义CPUE，其公式为

式中，CPUE为每年的单位捕捞努力量，单位为t/（船·a）；Catch的单位为 t/a。

利用广义线性回归模型（Generalized Linear Model，GLM）对名义CPUE进行标准化。由于CPUE受环境因子（SST）、渔场空间位置（渔场重心的经度、纬度）等因素影响，所以将其作为解释变量，CPUE作为响应变量，去除空间和环境因子对CPUE的影响，以名义CPUE和标准化的CPUE作为柔鱼资源丰度的指标。

渔场重心计算公式为

式中，Lon、Lat分别为渔场重心位置的经度和纬度；Ci为i渔区的产量；Loni、Lati分别为渔区i的中心经纬度位置；K为渔区的总个数。

GLM的表达形式为：glm[CPUE]～SST +Lon+Lat+ε。其中，SST为年平均海表面温度，ε为误差项。

2.3 模型的构建

利用GM（1, 1）模型对北太平洋柔鱼资源丰度进行灾变预测[15]。具体方法如下：

设原始序列X=（x(1)，x(2)，···，x(n)），取序列的第3分位数Q3（大于75%的值）为上限灾变值 ξ1，取序列第1分位数Q1（小于25%的值）为下限灾变值 ξ2，将序列X中大于或小于 ξ的点认为是异常值的点[16]，将其组成上限或下限灾变序列Xξ=（x[q(1)]，x[q(2)]，···，x[q(m)]），设灾变日期序列Q(0)=（q(1)，q(2)，···，q(m)），对其建立GM（1, 1）模型进行灾变预测。GM（1, 1）模型计算方法如下[17]：

令Q(1)=（x(1)(1)，x(1)(2)，···，x(1)(n)）为灾变日期序列Q(0)的一次累加生成序列。其中令Z(1)=（Z(1)(2)，Z(1)(3)，···，Z(1)(n)）为Q(1)的紧邻均值生成数列，其中，k=2,3,···,n，得到GM（1, 1）模型的基本形式为：Q(0)(k)+aZ(1)(k)=b，a为发展系数。

2.4 模型的检验

灰色灾变GM（1, 1）模型的验证主要包括残差检验与后验差检验[18]。残差为原始值与预测值之间的绝对误差与相对误差；后验差检验是按照残差的概率分布进行检验，记原始序列标准差绝对误差序列标准差方差比小误差概率

残差越小，模型精度越好。方差比C越小，说明S2小，S1大，即残差方差小，原始数据方差大，残差比较集中，摆动幅度小，原始数据比较分散，摆动幅度大，所以模型模拟效果更好。模型精度检验等级见表1。

发展系数a反映主因素与各子因素之间的协调程度，当-a≤0.3时，模型可用于中长期预测；当0.3＜-a≤0.5时可用于短期预测[19]。

表1 灰色系统理论模型精度检验[18]Table 1 Precision test accuracy of grey model[18]

3 研究结果

3.1 CPUE情况与灾变点

1995-2017年 各年的名义平均CPUE与标准化CPUE见图1。经过GLM模型标准化的CPUE与名义CPUE的变化趋势基本一致但波动较小。2007年名义 CPUE 值最大，为 4.197 t/（船·a），2009 年名义CPUE 最小，为 1.319 t/（船·a）；2017 年标准化的CPUE 值最大，为 2.928 t/（船·a），2002 年标准化的CPUE 最小，为 1.446 t/（船·a）。

根据灾变值的划分（下限灾变值为2.13 t/（船·a），上限灾变值为 2.39 t/（船·a）），以 GLM 模型标准化CPUE建立的灾变模型的上限灾变点为2005年、2006年、2007年、2008年、2010年、2017年，灾变序列为 11、12、13、14、16、23；下限灾变点为 1996年、1999年、2001年、2002年、2015年、2016年，灾变序列为 2、5、7、8、21、22。

3.2 模型的建立与检验

从模型精度来看，以GLM模型标准化CPUE建立的下限灾变预测模型的平均相对误差为15.32%，上限灾变预测模型的平均相对误差为8.19%（表2）。从发展系数a来看，下限灾变预测模型（a=-0.064）和上限灾变预测模型（a=0.117）均可用于中长期预测（表3）。从后验差检验结果（表3）来看，下限灾变预测模型（C=-0.342，P=1）和上限灾变预测模型（C=-0.341，P=1）的方差比C小于 0.35，小误差概率P大于0.95，模型精度等级为Ⅰ级。

根据表3的响应函数来看，下一次的下限灾变和上限灾变发生序号分别为34.86和25.20，在2016年以后的11.14 a（由34.86 a减去23.72 a得出），即2027年会发生资源丰度歉年；在2017年以后的4.0 a（由25.20 a减去21.20 a得出），即2021年会发生资源丰度丰年（表 3）。

4 讨论与分析

4.1 气候和环境变化对西北太平洋柔鱼丰度的影响

图1 1995-2017年西北太平洋柔鱼资源丰度变化Fig. 1 CPUE of Ommastrephes bartramii in the Northwest Pacific from 1995 to 2017

表2 灰色灾变GM（1, 1）模型预测值与真实值之间的相对误差Table 2 Relative error between the predicted value and the true value of the GM（1, 1）model

表3 灰色灾变GM（1, 1）模型的相关参数Table 3 Related parameters of the GM（1, 1）model

柔鱼作为短生命周期的种类，一生产卵一次，因此其渔获量的多少取决于当年的柔鱼补充量[20]。造成柔鱼资源丰度下降的原因是多方面的，如沈建华等[21]认为2001年低产的原因主要为：黑潮大蛇行、黑潮势力偏弱、北部水温过低、渔汛期流隔不明显和饵料不充足。陈峰等[22]认为2009年低产的原因主要为产卵场发生黑潮大弯曲，传统作业渔场有一个明显的冷水南下。唐峰华等[23]认为SST偏低导致的补充量减少，叶绿素a浓度分布与海流波动较大导致的饵料分布不稳定等现象也是2001年低产的原因之一。从大尺度气候变化来看，PDO是一种以10 a周期变化的太平洋气候变化，PDO指数定义为20°N以北北太平洋SST异常第一主分量。余为等[7]对不同气候模态下西北太平洋柔鱼渔场环境做出分析，认为PDO暖期会为仔幼鱼提供较好的摄食环境，易形成适宜的育肥场，导致CPUE较高，PDO冷期则相反。根据美国华盛顿大气与海洋联合研究所的数据（http://jisao.washington.edu/pdo/PDO.latest）结合灾变年份的分布来看（图2），基本上灾变年份符合PDO冷暖期的分布规律，即上限灾变年份位于PDO暖期，下限灾变年份位于PDO冷期。而且在冷暖期转换的年份（2001-2002年，2007-2008年，2014-2015年）灾变年份类型仍维持在转变期前的状态，说明PDO对柔鱼资源的影响具有一定的滞后性。此外，El Niño-La Niña事件是太平洋海域气候变率的最强信号，对柔鱼资源丰度也有较大的影响。如Chen等[6]认为La Niña事件发生时，柔鱼产卵场温度升高，不利于仔幼鱼的繁殖与孵化，导致资源补充量下降，反之El Niño事件时资源补充量上升。Yu等[24]通过构建栖息地模型发现La Niña事件下产生的暖水团，使适宜的叶绿素浓度和海表面高度距平范围扩大，导致栖息地面积扩张，而El Niño事件则相反。根据NOAA对El Niño和 La Niña事件的定义（http://ggweather.com/enso/oni.htm），即若Niño3.4区SST异常连续5个月滑动平均值超过+0.5℃认为发生一次El Niño事件；若连续5个月低于-0.5℃则认为发生一次La Niña事件。结合灾变年份分布来看（图3），基本上上限灾变年份发生在El Niño事件期间，下限灾变年份发生在La Niña事件期间。相比于 PDO 而言，El Niño 和 La Niña事件期间的灾变年份分布更为明显，如下限灾变年1996年对应La Niña事件，上限灾变年2003-2005年对应El Niño年，并且其滞后性比PDO更小。

4.2 灾变预测模型分析

灾变预测年份显示（表3），在4 a或者更短的时间内，北太平洋柔鱼的资源丰度依旧很高，这也符合今后持续的PDO暖期与El Niño事件的气候变化背景，但在 PDO 冷暖期与 El Niño-La Niña事件转化的年份，易发生资源丰度歉年，并在11 a左右的下一个周期开始仍处于低资源丰度的状态。所以建议北太平洋渔业生产企业和管理部门在PDO暖期与El Niño事件期间加大生产与捕捞努力量的投入，争取经济效益的最大化；相反，在PDO冷期与La Niña事件期间及时做好低产的准备，减少船只、人力、管理等方面由于多余投入而带来的经济损失。

图2 1995-2017年太平洋年代际涛动（PDO）指数时间序列变化及灾变年份分布Fig. 2 Time series of PDO and catastrophe year distribution from 1995 to 2017

图3 1995-2017年El Niño和La Niña事件变化及灾变年份分布Fig. 3 El Niño and La Niña event changes and catastrophic years distribution from 1995 to 2017

灰色灾变预测实质上异常值预测，而异常值的设置往往是根据人们凭经验确定的，具有很强的主观性，因此灾变阈值的确定十分重要。本研究采用统计学上的四分位数法确定灾变阈值的上下限，从模型的精度来看，灾变GM（1, 1）模型可以有效地预测灾变发生时间。从预测模型的结果看（表2），对于短期的预测值精度都较高，但随着预测时序的增长，会出现误差较大的点，证明灰色灾变预测模型仅适合短期的灾变点预测。此外，灾变预测不能预测连续发生的灾变点，如本研究的连续灾变年份包括上限灾变点2005-2008年，下限灾变点2000-2001年，2015-2016年。由于GM（1, 1）模型过于依赖样本序列，所以样本长度不同，序列起始划分不同，预测值会具有一定的不确定性，如样本太少会导致结果失真，样本太多难以达到所要求序列满足的光滑度。建议以后的研究采用多样本序列组合形式建立多种模型，考虑渔获量或捕捞努力量等指标，筛选出最优模型用于灾变长期预报。