刘振宇,孙伟富,陈 磊
(1.中南民族大学 资源与环境学院,湖北 武汉 430074;2.自然资源部第一海洋研究所 海洋物理与遥感研究室,山东 青岛 266061;3.中国科学院海洋研究所,山东 青岛 266071;4.中国石油大学(华东),山东 青岛 266580;5.中国科学院海洋大科学研究中心,山东 青岛 266071)
海表面盐度(Sea Surface Salinity,SSS)是表征河口冲淡水输入的关键指标,影响着海水密度、环流模式和物理过程[1-3]。海表面盐度的短期变化,会迅速改变海洋生物化学反应,进而影响海洋生态系统[4-5];而其长期变化则往往与气候变化存在关联[6-7]。因而,获取海表面盐度对海洋物理、生态环境以及大陆气候变化研究具有重要意义。
近年来,海表面盐度遥感技术因其广泛的空间覆盖和实时观测的优势,常用于替代现场观测。其中,微波辐射计SMOS(Soil Moisture and Ocean Salinity),Aquarius/SACD 和 SMAP(Soil Moisture Active Passive)被广泛用于获取全球开阔大洋的表面盐度[8-9]。然而,这些微波盐度卫星往往具备较低的空间分辨率(25~100 km)和时间分辨率(≥3 d)[10-11],加之L 波段微波辐射计受近岸射频(Radio Frequency Interference,RFI)干扰[12],无法提供渤海等近岸水体的海表面盐度数据。
近岸海域的有色溶解有机物(Coloured Dissolved Organic Material,CDOM)主要来源于陆地输入,因而它与海表面盐度通常存在线性负相关关系[5,13-15]。而CDOM 是光学遥感的可获取主要参数之一,因而光学遥感可获取存在这种负相关的近岸水体的海表面盐度。多波段组合线性回归,常被用于建立MERIS、MODIS、GOCI 等传感器的遥感反射率与海表面盐度之间的经验模型,以获取中国近海的海表面盐度[3,5,16-18]。近期研究发现,通过机器学习和神经网络建立模型,能够提高的盐度反演模型的精度[18-21]。
鉴于机器学习和神经网络在海表面盐度反演建模方面的优势,基于多层人工神经网络,本文建立MODIS-Aqua 二级产品提供的渤海遥感反射率与海表面盐度的经验模型,并评估模型精度。在此基础上研究渤海海表面盐度的时空变化。
渤海(北纬 37°07′~41°01′,东经 117°35′~121°10′)是中国的一个半封闭内海,面积约7.7×104 km2,平均水深18 m。它通过渤海海峡与黄海相连,主要由北部的辽东湾、南部的莱州湾、西部的渤海湾、中部的浅盆地和东部的渤海海峡组成,如图1所示。流入中国渤海的河流主要包括:黄河、海河、滦河和辽河,其中黄河径流量最大,约占总径流量的75%,导致渤海盐度的季节性变化[5]。
图1 研究区域和采样点分布Fig.1 Distribution of study area and sampling sites
共收集整理409 组由CTD(measuring conductivity,temperature,depth-based pressure)测量的海表面盐度数据,这些盐度数据分别采集自2007年1-9月、2009年5-10月以及2015年5月。根据《中国海洋调查规范》(GB/T 12763.3—2007),水深0.5 m 处CTD 测量盐度视作采样点的海表面盐度。采集样本的盐度直方图见图2(a),范围为26.64~33.03 psu,均值为30.98±0.86 psu,最低盐度的采样点位于黄河口附近,最高盐度则位于渤海湾。
图2 实测409 组盐度数据的直方图和匹配的MODIS-Aqua 遥感反射率Fig.2 Histogram of measured 409 sets of salinity data and matched MODIS-Aqua remotely sensed reflectance
MODIS-Aqua(以下简称 MODISA)卫星于2002年5月发射进入太阳同步轨道,自2002年6月运行至今,每隔1~2 d 覆盖地球表面。MODIS上设置了36 个波段,光谱范围400~1 440 nm,空间分辨率为0.25 km、0.5 km、1 km。NASA 提供的MODISA 二级产品(http://ladsweb.nascom.nasa.gov/data),包含空间分辨率为1 km 的412 nm、443 nm、469 nm、488 nm、531 nm、547 nm、555 nm、645 nm、667 nm、678 nm 等10 个波段的遥感反射率,鉴于667 nm 和678 nm 的大气校正精度相对较低,使用前8 个波段的遥感反射率。去除云和处理精度较低的像素(flags=4,5,6,9,10,15,20,21,23,30,http://oceancolor.gsfc.nasa.gov/ VALIDATION/flags.html);然后获取与实测数据时间相差不超过48 h[7]、空间位置匹配的MODIS 像素;以该像素为中心,统统计3×3 窗口平均,如果窗口内有效值超过5个,将窗口平均值作为实测数据的匹配结果,最终匹配得到164 组遥感反射率,见图2(b)。
1.3.1 渤海海表面盐度反演模型
由于遥感反射率与海表面盐度之间不存在直接关联,而是通过CDOM 建立起间接的关联,因而很难找到稳定的敏感于海表面盐度变化的特征波段和波段组合。因此,多波段线性回归模型是常用于反演渤海海表面盐度,针对MODIS 数据QING等[3]发展的模型主要使用3 个水色波段490 nm、560 nm、665 nm:
公式(1)中给出的波段组合还被扩展到GOCI数据[5,17]以发展渤海海表面盐度反演模型。YU等[16]认为MODIS 数据531 nm 和555 nm 波段遥感反射率比值是反演渤海CDOM 的特征波段,因而可用于获取渤海海表面盐度:
1.3.2 多层神经网络
多层神经网络(以下简称“神经网络”)是一种前馈式人工神经网络模型,用于模式分类、识别、预测和函数逼近。研究表明,具有一个或多个隐藏层和非线性激活函数的神经网络可以逼近非线性函数[22-23]。因而它适合用于建立存在间接关联的遥感反射率与海表面盐度之间的经验模型。对于神经网络建模而言,核心问题是选择最优的隐层数和神经元数。这通常受很多因素的影响,如输入和输出层的大小、训练样本的数量、待逼近函数的复杂性、激活函数的类型以及训练算法等等。鉴于遥感反射率与海表面盐度之间的间接关联,使用1~2 个隐层配置的神经网络。
使用MATLAB 提供的神经网络工具箱,随机选70%的数据作为训练数据集,剩下30%用于模型验证、测试;根据训练数据集中盐度的最大和最小值,通过线性拉伸对海表面盐度进行归一化。通过比对不同配置下神经网络的表现,最终确定模型配置为1 个隐层、30 个神经元,sigmoid 激活函数。
本文研究流程图3所示,利用1.2 节提到的方法用实测海表面盐度数据去匹配MODISA 二级产品,从而得到164 组MODISA 遥感反射率(412 nm、443 nm、469 nm、488 nm、531 nm、547 nm、555 nm、645 nm),然后分别采用多波段回归和神经网络建模,验证、比对模型精度,最后分析渤海海表面盐度的时空变化。
图3 本文研究的技术流程Fig.3 Technical flow of our research
利用本文的164 组匹配数据,采用多波段线性模型(公式(1)和(2)),建立渤海海表面盐度反演模型(公式(3)和(4))。这2 个模型与本文的神经网络的精度比对见图4。
图4 神经网络与多波段线性模型建模精度比对Fig.4 Accuracy comparison between neural network model and multi-band linear model
注意:神经网络模型只使用随机挑选的114 组数据(约占总数据的70%),其(图4(a))决定系数R2和均方根误差RMSE 分别是0.66 和0.39,均优于多波段线性模型的精度(0.39 和0.60)。事实上,海表面盐度与MODISA 各波段遥感反射率的统计相关性,受数据采集时间和采样区域影响,存在较大变化[24]。也就是说,很难找到一个敏感于盐度变化的固定波段组合。由式(3)和式(4)波段组合建立的模型,当盐度小于31 psu 时,拟合的散点出现较大偏差。这是导致其精度降低的直接原因。
利用空间分辨率为4 km 的MODISA 季度平均遥感反射率(三级产品),经本文的神经网络模型计算,得到2022年3-6月、6-9月、9-12月、2022年12月-2023年3月等,共4 个季度的渤海海表面盐度分布(见图5)。从空间分布来看,渤海湾和莱州湾的海表面盐度较低,渤海海峡的盐度相对较高,这与之前的研究结果一致[3,5,16-18]。这4个季度渤海的平均盐度分别为 31.25±0.40 psu,31.52±0.39 psu,31.13±0.60 psu、30.89±0.53 psu。黄河流域的雨季主要集中6-9月,然而这个季度渤海的海表面盐度最高,表明黄河流域的降水对黄河冲淡水影响存在一定时间的滞后。图5(c)的指示箭头可见,黄河冲淡水对莱州湾和渤海湾盐度的影响,同时也会随着渤海冷流南下致使山东半岛沿岸的盐度降低;至第4 季度(图5(d)),南下的黄河冲淡水向渤海海峡和北黄海方向扩展,并持续影响莱州湾,导致莱州湾以及整个渤海盐度降到最低。
图5 2022年3月-2023年3月渤海海表面盐度季度分布Fig.5 Quarterly distribution of sea surface salinity in the Bohai Sea from March 2022 to March 2023
过拟合是指模型学习了训练集的噪声或细微特征,导致模型泛化能力下降,当应用于测试数据集时精度显著下降;样本数量少、模型复杂度高(隐层的数和神经元数)、训练时间过长等是导致网络过拟合的主要原因[25]。为避免过拟合,本研究对海表面盐度进行归一化处理,同时通过比对训练数据集(占数据集70%)和测试数据集的损失函数变化。当训练数据集损失函数降低,而测试数据集损失函数上升,表明模型可能已过拟合。
鉴于遥感反射率与海表面盐度通过CDOM 存在间接关联,首先考虑2 个隐层的神经网络。相比1 个隐层,2 个隐层(神经元数目分别为15 和25,(记做[15,25])的网络能提高模型精度(R2≥0.70),应用于MODISA 的遥感反射率数据时,得到结果与前文给出模型比对见图6。渤海中部和山东半岛沿岸展现出高盐度,甚至超过渤海海峡的盐度,这与前期研究[3,5,16-17]结果矛盾,而采用单隐层的网络模型结果与先前研究一致。这表明2 个隐层网络已过拟合。因而文中选用单隐层神经网络。
图6 2022年渤海平均海表面盐度分布Fig.6 Distribution of mean sea surface salinity in the Bohai Sea in 2022
由于海表面盐度与CDOM 呈负相关,总有机吸收系数(adg=CDOM 吸收系数ad和有机碎屑的吸收系数ag,单位m-1),也被证实与海表面盐度存在良好负相关[14,26]。此外,作为海水中有机颗粒物浓度表征的颗粒有机碳(Particulate Organic Carbon,POC),在近岸区域也像CDOM一样受陆源输入影响。因此,利用MODIS 二级产品中443 nm 总有机吸收系数adg(443)和颗粒有机碳,与实测海表面盐度数据匹配后,并剔除异常数据(adg(443)>1 m-1;POC>600 mg/m3),通过线性回归建立两者之间的经验关联,结果见图7。
图7 海表面盐度与adg(443)和POC 的线性回归模型Fig.7 Linear regression model of sea surface salinity with adg(443)and Poc
adg(443)和POC 与海表面盐度的R2分别为0.24 和0.33,虽低于神经网络模型的0.66,但也表明MODISA 二级产品中的这2 种数据可作为获取海表面盐度辅助数据。至于这2 种产品的自身精度,以及它对图7 线性回归的影响,不属于本文的研究内容。
用164 组实测海表面盐度和匹配的MODISAqua 遥感反射率产品,通过训练神经网络建立了渤海海表面盐度的反演模型,其R2和RMSE 分别为0.66 和0.39,均优于多波段线性模型(R2=0.39,RMSE=0.60)。利用该模型,绘制了2022年3月-2023年3月4 个季度的海表面盐度分布,占渤海径流总量75%的黄河冲淡水是影响渤海海表面盐度的重要因素;3-9月,黄河冲淡水主要扩散至莱州湾和渤海湾,导致其盐度降低;9-12月,除持续影响渤海湾和莱州湾,黄河冲淡水随着渤海冷流沿山东半岛南下,致使沿岸以及整个渤海的盐度降低,随后向渤海海峡和北黄海方向扩散。以上结论为渤海海表面盐度和黄河冲淡水研究提供技术支持。