渤海是北半球纬度最低的结冰海域, 典型的一年冰海区, 每年11月末到次年3月初总会出现不同程度的结冰现象(王萌等, 2016), 受海陆环境影响, 有很强的季节性变化特征(Jones et al, 2016)。海冰暴发的时间、空间差异以及海冰强度极端现象使渤海成为受海冰灾害影响的主要区域之一, 例如在2010年为重冰年, 海冰几乎覆盖整个渤海, 直接经济损失达63.18亿元(孙劭等, 2011), 为避免或减小海冰灾害就必须掌握海冰发展的时空变化规律; 同时海冰又是一种重要的淡水资源, 成为缓解环渤海区域淡水资源危机的新方式、新途径(史培军等, 2002), 海水淡化的实现同样需要确切掌握渤海范围内海冰时、空分布区域差异。辽东湾是渤海海冰灾害最严重的海域, 精确预测海冰时空分布规律, 对于预防海冰灾害和应用具有重要的现实意义。

渤海海冰一直是国内学者研究的热点。在海冰监测方面, 赵进平等(2003)根据边缘区的温度及盐度特征可以将北冰洋海冰边缘区分为暖水边缘区和冷水边缘区, 为渤海海冰研究提供借鉴; 吴龙涛等(2006)根据海冰与海水光谱特征差异, 实现利用遥感数据对海冰密集度的提取; 郭衍游等(2010)利用海冰表面温度反演算法进一步区分海冰和泥沙等悬浮物。对于海冰规律研究, 张云吉等(2007)认为渤海海冰是对华北地区中短期气候变化的响应, 从20世纪50年代到90年代冰情等级减弱; 李彦青等(2013)利用海冰离岸距离建立侯平均时间序列, 在小时间尺度上分析2001—2011年渤海海冰的变化。至于影响海冰发展的因素, 一些学者认为日均气温积温与渤海逐日冰面积存在较好的相关性, 是影响海冰发展的主要因素(李春花等, 2008; 武晋雯等, 2009); 顾卫等(2002)得出–100℃等值线与渤海常年冰情分布范围大致相当; 唐茂宁等(2015)指出西伯利亚高压和西太平洋副高也对渤海海冰强度产生影响; 此外部分学者认为渤海海冰发展还与寒潮(朱晨玉等, 2014)、厄尔尼诺(白珊等, 2001)等有一定相关性。海冰预测一直是海冰研究的难点, 建立海洋环境数据库, 模拟分析大气、海洋之间的关系, 是海冰预测的基础(李春华等, 2003); 余加艾等(1999)基于溢油与海水、海冰、大气的相互作用原理, 建立了渤海结冰海区溢油行为数值模拟; 张娜(2012)在考虑水动力、海温、盐度、海冰热力学等基础上建立非结构化网格FVCOM三维海洋模型对实现对海冰中长期数值模拟, 数据所需要素过多一定程度上增加模型使用的难度; 周须文等(2015)采用主成分回归分析(PCR)法, 建立逐旬滚动的渤海冰情的定量预测模型, 但在预测的准确率及稳定性上存在缺陷。

以往学者对海冰的研究多集中于对海冰提取方法、年际变化规律, 影响因子的研究, 对海冰的预报、动态生消变化研究相对较少, 仍需要进一步加强。辽东湾位于渤海北部, 东西北三面被辽宁省环绕(图 1), 受高纬度影响, 冰情重于莱州湾和渤海湾, 是渤海冰情最严重的海区, 周围有辽河、大凌河等水系注入, 盐度较低, 泥沙淤积, 海水深度由湾口向渤海内部随离岸距离逐渐加深, 最大深度达到32m, 受温带季风气候影响, 冬季寒冷, 每年冬季都会出现不同程度的结冰现象(李秉天等, 2015; 马伟伟等, 2016)。环渤海地区又是我国主要的经济开发区之一, 经济发达, 航运繁忙, 资源丰富, 海水养殖兴盛, 致使海冰灾害成为渤海最严重的海洋灾害之一。本文以MODIS为数据源, 辽东湾为研究区, 利用2010/2011年度到2014/2015年度海冰结冰情况, 计算出辽东湾海冰结冰概率, 建立营口气象站积温与结冰概率的关系, 对辽东湾2015/2016年度海冰发展及空间分布进行模拟, 提出一种局部条件下海冰预测的新思路, 为海冰灾害预防, 淡水利用, 精确预报海冰发展的空间分布奠定基础。

1 数据来源

本文研究使用的数据主要有两类, 海冰结冰信息数据和气象数据。其中海冰结冰信息是利用EOS-AQUA和EOS-TERRA影像提取, 遥感影像来源于美国国家航空航天局(网址: https://ladsweb.nascom.nasa.gov/search/), MODIS数据可对同一地区一天上午、下午2次重复观测, 空间分辨率为250m、500m、1000m之间, 对海冰与海水的判别精度较高, 成为国内学者对渤海海冰监测的主要数据之一, 与常用的Landsat系列卫星相比在时间分辨率, 及覆盖范围上有明显优势(郭衍游等, 2014; 奚萌等, 2017)。

本文选取2011—2016年无云或少云(遥感影像局部区域被云覆盖但不遮挡研究区范围; 或被薄云覆盖, 透过云雾能够判断研究区内海冰冰情时仍视为有效数据)的MODIS 1B产品, 共177景(可用数据), 年均35景, 涵盖初冰期到盛冰期(60天左右)全过程, 可用率在30%—60%, 空间分辨率为500m, 其中红色和近红外波段分辨率达到250m。

气象数据主要是2011—2016年, 环渤海气象站(营口站)日均值数据, 来源于中国气象科学数据网站(网址: http://data.cma.cn/)中国地面国际交换站气候资料日值数据集。

2 研究方法 2.1 海冰提取方法

对MODIS遥感影像进行大气校正、几何校正、裁剪等预处理, 利用海冰在可见光波段(0.38—0.76µm)反射率较高, 短波红外波段(1.57—1.65µm)反射率较低(潘淑杰, 2008)的反射率特征, 使用近红外(R)、红光(G)、绿光(B) 3个波段假彩色合成, 在假彩色影像上海冰呈现亮蓝色, 海水呈暗蓝色, 陆地偏红色的特征, 色差鲜明, 利用目视解译方法, 基于ENVI和ArcGIS软件对海冰与海水分界线矢量化, 进而利用海冰边缘线对渤海边界裁剪, 实现对海冰的定量提取(张启文, 1986)。

2.2 结冰概率

辽东湾海冰结冰程度受温度、盐度、风力、海水深度及周边地形等多种因素综合影响, 且不同影响因子对海冰影响程度不同, 若精确量化海冰每一种影响因子(如图 2a)对海冰发展程度进行预测具有一定难度。辽东湾是内陆海湾, 三面环陆, 海陆环境单一, 区域环境内气候稳定, 海冰发展存在一定规律, 随气温的降低, 辽东湾海冰由湾底向湾口, 由两岸向中轴线方向发展。在此基础之上, 可以认为在海冰冰情等级相似的年份, 各种影响因子的综合作用产生效果相似, 可以忽略各种海冰影响因子的差异, 直接观察海冰发展变化(图 2b), 利用已知年份海冰发展规律模拟预测未知年份海冰发展变化。为进一步量化海冰发展规律, 本文提出结冰概率的概念, 即区域范围内每一点的结冰天数与区域总结冰天数的比值, 利用区域内结冰概率的高低来展示海冰发展变化。

2.3 结冰概率计算方法

为进一步计算辽东湾内部结冰概率的具体分布, 首先需要获取研究区内结冰信息, 本文将2010/2011年度—2014/2015年度海冰日结冰数据统一排序(1、2、3……k期), 将研究区划分为100×100的格点, 格点水平分辨率为4.5km(即每间隔9个像元有一个格点), 格点在海冰范围内时, 格点属性赋值为1(判断为海冰), 否则为0, 计算公式如下,

式中, i、j、k分别是格点的行号(i=1, 2, …, 100)、列号(j=1, 2…100)及海冰期数, Pt_ijk是第i行, j列, k期格点的属性值, U_k表示k期内海冰结冰区域。

式中, Pt_ij为第i行, j列格点结冰的总期数, 即为2010/2011年度—2014/2015年度间每个格点的日结冰次数的属性值, 将格点属性值累加, 根据格点属性值利用克里金差值等方法分析海冰分布特征。

进一步研究渤海海冰在空间分布规律, 本文根据2010/2011年度—2014/2015年度的结冰信息, 计算出辽东湾不同海域的结冰概率, 计算公式如下(式3)。

式中, Ps_ij是第i行, j列单个格点的结冰概率, Max(Pt_ij)是研究区100×100格点中最多的结冰期数, 每年海冰由高结冰概率区域按一定的轨迹向低结冰区域发展, 结冰概率越高表明该地易于结冰, 冰封日期越长。利用克里金插值方法对概率点插值, 得到研究区结冰概率分布(图 3), 为下一步利用积温与结冰概率的关系动态预测海冰发展奠定基础。

2.4 营口积温与结冰概率

在众多海冰影响因子中气候因素起主要作用(张启文, 1986), 特别是气温的变化直接影响海冰的发展(赵春雷等, 2014)。营口站积温与辽东湾海冰的发展过程有较好的相关关系, 每年从12月初到次年2月, 随营口积温逐渐降低, 辽东湾结冰范围由北向南、由沿岸向中轴线方向逐渐增加。本文以气温作为影响辽东湾海冰发展主导因素, 将营口冬季12月1日起积温数据与辽东湾结冰概率相结合, 提取2010—2015年海冰日外缘线及其积温条件下的结冰概率, 相关关系如下(图 4)。

从总体分析, 结冰概率与营口积温有较好的相关性, 随营口积温的减小, 海冰结冰范围由北向南, 辽东湾沿岸向海湾中轴线方向, 由高结冰概率区域向低结冰概率区域发展。绿色曲线是概率点的内插线, 反映结冰规律随营口站积温的起伏变化及总体趋势; 经过概率点非线性拟合的曲线(红色实线)过高的估计了结冰的概率, 特别是在纬度相对较低区域误差较大, 在低结冰概率的区域出现结误判现象, 当完全按照最低结冰概率构建模型时, 预测概率又过低的现象; 为更加符合结冰实情, 本文以同积温条件下以较低概率为基准, 在内插线下折点与同积温条件下结冰概率中点处及结冰概率的始终点, 选取25个采样点(图 4)进行曲线拟合。由表 1可知, 从结冰概率与及稳定整体关系来看, 指数拟合曲线R²=0.789, 拟合优度高于线性拟合及二次型拟合曲线, 基于此对采样点使用指数函数模型拟合(蓝色实线), 得到R²=0.943, 拟合优度较高, 拟合公式如下,

式中, P_pre是估计的结冰时最低达到的最低结冰概率, T为营口站冬季(12月1日起)积温(℃), 当Ps > P_pre的区域为结冰区。

3 结果与分析 3.1 气候影响因子

营口站积温年际变化较大(图 5), 在初冰期(12月份)积温差距较小, 之后随温度的降低积温差距逐渐增大, 2010/2011—2015/2016年度最低积温出现在2012/2013年度为–732.7℃, 2013/2014年度与其它年度积温相比较高, 最低积温在–432.9℃, 其中2010/2011和2011/2012年度积温变化相似。温度是影响海冰的主要因素, 营口站积温与辽东湾海冰的发展存在较高的相关性(王萌等, 2016; 白珊等, 2001), 这使得利用营口站积温预测海冰的发展及空间变化成为可能。

3.2 海冰预测

初冰期至盛冰期阶段是海冰产生到致灾的主要时期, 这个时期内海冰发展迅速, 逐渐蔓延到整个辽东湾, 成为海冰预监测、预报的关键时期。本文以2010/2011年度—2014/2015年度结冰信息为基础, 对2015/2016年度辽东湾初冰期至盛冰期阶段进行预测, 构建结冰概率与营口气温的关系模型, 探讨海冰局部预测的可能性。

图 6由近红外、红光、绿光波段假彩色合成, 土壤呈偏红色, 海冰呈亮蓝色, 海水暗蓝色的特征, 海冰和海水易于区分。2015/2016年度海冰从2015年12月4日开始结冰, 到2016年2月2日结冰面积达到最大, 2016年3月8日完全消融, 从总体上看, 海冰预测结果与海冰实际情况较为吻合, 预测海冰外缘线与海冰假彩色影像结果相近(图 6)。图 6a是海冰发展的初期, 海冰主要由辽东湾湾底向湾口方向发展, 海冰发展与纬度平行, 图 6b海冰由湾底向湾口及海湾两侧向中轴线两个方向发展, 图 6c海冰面积达到2015/2016年度辽东湾最大值, 海湾基本被海冰覆盖。

3.3 精度验证

进一步验证海冰提取结果, 将2016年1月10日、1月21日和2月2日MODIS影像目视解译所得的分类结果作为真实海冰数据, 概率关系模型提取的海冰分类结果设定为预测值, 利用混淆矩阵方法分析海冰提取精度, 具体结果如下(表 2、表 3、表 4)。3个时期的分类总体精度分别为95.0697%、93.0731%和93.6683%, Kappa系数K_hat分别为0.9142、0.8862和0.8957, 总体精度在93%以上, 说明分类结果有较高的有效性, 利用营口站积温与海冰概率的关系可以预测辽东湾海冰的动态变化, 为其它地区小区域环境内海冰预测提供参考。

海冰预测结果与真实分布情况存在一定差异(图 7), 且主要集中在海冰与海水的交界区域, 总体可以来看将海水误判为海冰的面积要大于海冰误判为海水的面积, 1月10日、1月21日及2月2日误判为海冰的面积分别占研究区的4.23%、5.3%和7.86%, 随海冰面积的增加误判区域越大。

为进一步检验海冰预测方法的稳定性, 利用相同的方法以2003/2004—2007/2008年度辽东湾的结冰信息为数据源计算出结冰概率信息, 基于营口站积温与结冰概率模型, 数值模拟辽东湾2008/2009年度结冰情况, 间隔10天左右选取6景加以展示(图 8)。总体来看, 2008年12月14日、12月31日、2009年1月2日、1月11日预测效果优于2月18日、2月3日, 预测海冰与实际冰情利用混淆矩阵验证精度分别为98.7973%、96.0773%、96.3523%和90.6708%, 1月18日、2月3日海冰预测结果与实际海冰分布差异较大。2009年辽东湾冰情等级属于偏轻冰年, 与往年同积温条件下结冰区域明显减小, 1月18日前后, 其它因素的影响大于温度对海冰的影响, 出现积温持续减小(温度降低)而海冰面积没有增加甚至减小的现象, 海冰增加和消融波动起伏较大, 增加了预测难度, 致使1月18日模拟效果较差。2008/2009年度海冰融冰期较往年提前, 2009年2月3日辽东湾已进入融冰期, 海湾西侧出现融冰现象, 本文海冰预测方法主要针对初冰期至盛冰期阶段的海冰发展过程, 因此与预测结果有所差异。基于2008/2009和2015/2016年度海冰模拟验证结果, 利用营口积温与结冰概率构建的模型, 一定程度上, 能够有效实现海冰动态预测。

4 讨论

积温是影响辽东湾海冰发展的首要因素, 营口站积温的高低与辽东湾海冰发展程度具有一定的相关性(李宁等, 2008; Zhang et al, 2016), 王萌等(2016)利用积温实现了对海冰面积的预测, 但海冰面积在辽东湾海区的具体分布难以确定, 成为海冰预报的难点。本文利用辽东湾海冰多年结冰实情计算得到的结冰概率, 可以直观的反映该区域海冰发展的轨迹、区域范围内海水结冰可能性高低, 通过研究辽东湾外缘线所处结冰概率的高低随营口站积温变化, 最终确定了辽东湾海冰随积温的动态分布, 一定程度上实现了辽东湾海冰的动态预报, 为其它小区域环境海冰预报提供方法参考。

基于海冰结冰概率与营口积温关系对海冰进行实时预测, 基本思想是利用近似条件下海冰发展轨迹预测未知轨迹, 预测精度受基础数据的选择影响, 即已知海冰发展年份与预测年份外环境越相近, 则预测精度越高。国家海洋局将冰情等级划分为重冰年、偏重冰年、常冰年、轻冰年、偏轻冰年五个等级, 在同一冰情等级年份中辽东湾积温等主要影响因子差异较小, 海冰发展轨迹相近, 可以利用同一冰情等级下结冰概率与营口站积温的关系实现对辽东湾海冰的预测。郑冬梅等(2015)学者海冰冰情等级具有2—7年的周期性变化, 海冰冰情等级与初冰期时间有一定关系(陈艳红等, 2015), 这些都为预测海冰冰情等级及基础数据的选择提供参考。模型的建立也存在一定局限性, 采样点选取采用人工随机选取, 其严密性存在不足, 会对建模精度产生一定影响, 大量实验及经验的基础上建议, 在概率点内插线下部的折点与同积温条件下最低概率点的中值附近、始点及终点处选取, 一定程度上可以避免预测值过高或过低的现象; 模型的构建是建立在数学统计基础之上, 使用非线性拟合方法, 在机理上的解释需要进一步探讨; 李冬玲等(2016)对西格陵兰Disko湾的研究表明, 沿岸洋流的冷、暖水团强度变化也是影响海冰的因素之一, 渤海湾沿岸洋流对结冰概率模型预测精度的影响尚未探究; 本文在对辽东湾海冰实时冰情预测, 主要是基于积温因素, 较为单一, 且在同一积温条件下辽东湾东、西两侧结冰概率存在差异, 在后续工作中将进一步对辽东湾进行区域划分, 引入其它海冰影响因子, 精确量化、完善不同地理环境下海冰预报方案。

5 结论

本文利用2010/2011年度到2014/2015年度营口站每年12月到次年2月的日气温气象数据对辽东湾地区海冰影响因子分析, 发现辽东湾积温差异较大, 成为影响海冰发展的主要因素; 根据2010—2015年辽东湾海冰结冰信息计算得出该海域结冰概率, 建立营口气象站积温与结冰概率的关系模型, 可以实现对辽东湾海冰的动态预测, 通过混淆矩阵对海冰提取结果检验, 总体精度在90%以上; 为了检验海冰预测方法的稳定性, 利用相同预测方法, 以2003/ 2004—2007/2008年度海冰为数据源对2008/2009年度海冰发展变化进行动态预测, 预测结果与实际冰情较为吻合, 为渤海及其它局部区域条件下初冰期至盛冰期海冰预报提供参考。