基于深度神经网络算法的水体透明度反演方法

喻臻钰,杨 昆,罗 毅,*,商春雪,赵 磊

1 云南师范大学地理学部,昆明 650500 2 云南师范大学西部资源环境地理信息技术教育部工程研究中心,昆明 650500 3 云南师范大学教务处,昆明 650500 4 云南师范大学信息学院,昆明 650500

湖泊既是重要的生态资源,同时也是一系列社会经济过程的载体,湖泊水生态环境不仅决定着区域生态环境质量,也是城市的可持续发展的重要自然资源保障[1- 3],湖泊水生态环境能够反映该地区环境状况[4- 6]。水体透明度(Secchi depth, SD)是指水体的澄清程度,水体透明度会随着水体中的悬浮物和胶体浓度增加而降低,在富营养湖泊中,以藻类为主的悬浮物越多,其透明度越低[7]。水体透明度通过微粒和溶解物质的光吸收影响上层水体的传热,水体透明度观测往往是湖泊水质监测的组成部分[7-8]。因此,水体透明度是衡量水质标准的重要组成部分,也是评价水体富营养化和水生态健康的重要指标[8],同时在调节水层和初级生产力方面发挥着重要作用[9]。近年来,由于外来物种入侵、污染物点源排放、不透水面增加、营养物质负荷过重以及由此造成的有害藻类大量繁殖等综合影响,滇池水体透明度出现了明显的波动[10-11]。水体透明度与水质关系密切[12-13],常用作度量湖泊总体水质的直观指标。

原位观测是掌握湖泊水体透明度常用的方法,具有数据时间尺度连续、精度较高等优点,但存在传感器布设和维护成本高、数据空间尺度离散等缺点[14- 16]。遥感技术的发展为区域环境监测提供了有效支撑,遥感影像的地表反射率数据为湖泊水体透明度的时空变化检测提供了基础[9]。Kloiber等[15]以Landsat 系列影像为数据,提出了Pearson相关系数矩阵和多元逐步回归的透明度反演经验模型,Olmanson等[17]以此模型对明尼苏达州的10000多个湖泊水体透明度进行了反演,并分析了其20年间的时空变化特征,表明Landsat系列影像数据在透明度反演中具有较高精度；Wu等[9]以Landsat和MODIS为数据,以多元线性回归模型对鄱阳湖的透明度反演结果表明,MODIS效果优于Landsat；马建行等[18]及Knight等[19]以HJ- 1卫星CCD和MODIS为数据,基于灰色关联度和多元回归反演模型,对吉林省中西部湖泊透明度进行反演和对比,结果表明MODIS反演精度高于HJ1A-CCD。MODIS影像虽然空间分辨率不高,但在水体透明度反演中表现较好,且其具有较高时间分辨率,为较大尺度的湖泊水体透明度实时变化监测提供了可能。

在水体透明度反演模型中,半分析方法和经验方法较为常用,其中半分析方法考虑了透明度与水体光学特性之间的关系,具有较高反演精度,而经验方法对云污染、大气影响不敏感,具有更强鲁棒性[18]。随着机器学习领域的发展,深度学习方法的提出,遥感影像反演精度得以大幅提升,现已被用于地表温度反演[20]、海表温度反演[21]、叶绿素a浓度反演[22]等,但尚未检索到有关深度学习理论应用于水体透明度反演的报道。

地表反射率与水体透明度的数据均具有混沌特征,多项式拟合的方法并不能最大限度的挖掘其有效特征。为此,在前期相关工作的基础上[11,23- 29],以原位监测和MODIS遥感影像为数据,基于深度神经网络方法,对滇池2001年1月1日—2018年12月31日的水体透明度进行反演,并利用地理空间分析方法对其时空变化特征进行探究,填补了深度学习在水体透明度反演的空白,有效提高了常用水体透明度反演模型的精度,为开展湖泊水质研究提供了新的思路和方法。

1 研究区与数据来源

1.1 研究区概述

滇池位于云贵高原中部,呈南北向分布,是全国第六大淡水湖,湖面海拔约1886 m,面积约330 km2,平均水深约为5 m,属于半封闭型湖泊,仅有西南部的海口为出水口,处于亚热带高原西南季风气候区,气候变化主要受西南季风和热带大陆气团交替控制。同时,滇池是昆明地区灌溉、调畜、受纳的主要水体,也是城市发展的基础条件,对滇池流域的气候调节具有重要作用。研究发现,近30年来,在城镇化的快速推进和全球气候变暖双重胁迫下,以湖泊表面水温上升、蓝藻水华频繁、富营养化加剧为代表的水质恶化状况明显加剧[11,23- 24,26,28],湖泊水环境不容乐观。另外,随着国家“一带一路”战略的提出,云南省正深度融入国家战略,作为云南省省会的昆明市,已经成为面向南亚、东南亚重要的辐射中心,控制与改善滇池水生态环境是落实和推进国家“一带一路”战略的环境基础。为此,选择滇池作为研究区具有区域特色(高原城市型湖泊)与时代特色(“一带一路”战略),研究区及监测站点地理位置分布如图1所示。

图1 研究区位图Fig.1 Study area

1.2 数据来源

时间区间为2001年1月1日—2018年12月31日,主要数据包括水体透明度实际观测数据和MOD09GA遥感影像。其中,实际观测值为2001—2010年10个站点(白鱼口、草海中心、滇池南、断桥、观音山东、观音山西、观音山中、海口西、晖湾中、罗家营)的日监测值和2018年2个区域(草海、外海)的月监测值,来源于云南省环境科学研究院；MODIS数据来源于NASA LAADS Web (http://ladsweb.nascom.nasa.gov),时间分辨率为每日,空间分辨率为500 m,利用MRT (https://lpdaac.usgs.gov/tools/modis_reprojection_tool)重采样为500 m分辨率影像,并将其重投影、镶嵌并替换有云影响的像元值为Null,遥感影像缺失数据的时间分布如图2所示,对于红色标注的缺失像元以线性插值法补齐。

图2 MODIS遥感影像缺失数据的时间分布Fig.2 Time distribution of missing value in MODIS remote sensing images

2 实验方法

2.1 灰色关联度分析(Grey Relational Analysis, GRA)

灰色系统理论利用关联度度量因素之间的关联性,不同灰色关联度方法的计算存在一定局限性[30-31]。本文采用一般关联度、绝对关联度、斜率关联度三种方法得到灰色关联度的相关系数,考虑到三种关联度的计算无法体现因素的负相关特性,因此采用Pearson相关系数分析负相关特性。

2.2 长短期记忆神经网络(LSTM)估算模型

每个Cell包括输入门、遗忘门和输出门(如图3所示),每个信息进入LSTM的网络,可以根据规则判断是否有用,只有符合规则的信息才会留下,不符的信息则通过遗忘门被遗忘,这对于具有长期序列依赖问题的数据非常有效[32]。

it=σ(Wxixt+Whiht-1+bi)

(1)

ft=σ(Wxfxt+Whfht-1+bf)

(2)

ot=σ(Wxoxt+Whoht-1+bo)

(3)

ct=ft⊙ct-1+it⊙tanh(Wxcxt+Whcht-1+bc)

(4)

ht=ot⊙tanh(ct)

(5)

其中：⊙指矩阵逐元素点乘；bγ是各层输出的偏差向量,例如bi是输入门限层的偏差向量,bf是遗忘门限层的偏差向量；σ(x)是激活函数；Wαβ是对应层的权重矩阵,例如Wxf是输入层到遗忘门限层的权重矩阵,Whi是隐藏层到输入门限层的权重矩阵,Who是隐藏层到输出门限层的权重矩阵；ct用来更新细胞状态。由式(4)可知,遗忘门ft控制有多少上一时刻的记忆细胞中的信息ct-1可以传输到当前时刻的记忆细胞中；输入门it控制有多少信息可以流入记忆细胞ct中；而输出门ot控制有多少当前时刻的记忆细胞ct中的信息可以流入当前隐藏层ht中[33]。

2.3 湖泊水体透明度反演模型(COD-LSTM)

数据包括MOD09GA和透明度的实测数据,其中MOD09GA为2001—2018年的日数据,实测数据为2001—2010年10个站点的日监测数据和2018年2个区域的月监测数据。以2001—2010年的实测数据与对应时间的MOD09GA数据进行建模,并以此模型对滇池2001—2018年的水体透明度进行估算,以2018年的实测数据对估算数据进行测试,技术路线如图3所示。

(1)相关性分析(Correlation Analysis, C)。以灰色关联度对MOD09GA的7个波段组合进行分析,得到排名前7的波段,即b5/b6、1/b4、b3/b6、b2/b6、b1/b6、b2/b3、1/b5,并组合为数据集X,其中bi为第i波段。

(2)异常值处理(Outlier Processing, O)。根据数据集的特征,临近时间点的数据存在较为紧密的线性关系,所以使用线性插值法对缺失值进行插补。

(6)

异常值以3σ准则进行检验,并将大于μ+3σ的离群值以99%分位数插补,小于μ-3σ的离群值以1%分位数插补,其中μ为均值,σ为标准差,如式(6)所示,此插补方法即盖帽法。

(3)去噪处理(Denoising Processing, D)。实际观测值与遥感影像均不能排除仪器设备及人为因素的干扰,本文使用高斯低通滤波器对数据进行去噪处理。高斯低通滤波器(Gaussian Low Pass Filter)是一类传递函数为高斯函数的线性平滑滤波器,其传递函数如式(7)。

地理定位信息系统的主要原理是将与地理相关的内容进行电子化的转换，使其以动态的图形的形式呈现出来。这种信息系统在旅游管理中的应用主要体现在对游客所在的位置进行实时监控[1]。旅游管理人员通过这样的方式掌握游客的位置信息，可以为游客提供安全保障。所以，在旅游管理中应用地理定位信息系统，不仅能使游客感受到更加优质和完善的服务，还能为促进旅游管理信息化建设提供推动力。

(7)

其中δ为标准差,x为待滤波数据,G(x)为滤波后的数据。

(4)估算模型构建。经过以上步骤处理的数据作为LSTM模型的输入,数据总数为36520行(10×365×10+2×1×10),其中80%(29216行)作为训练集(DataT),20%(7304行)作为验证集(DataV)。模型参数设置：经过以上3个步骤的7个波段组合数据X作为模型的输入(输入维度为7),隐含层为10,学习率为0.0006,每批次大小(Batch size)为7,待估算的水体透明度数据为模型的输出(输出维度为1)。

(5)结果与分析。对输出的水体透明度以Theil-Sen估算进行时序变化特征分析,以反距离加权(Inverse Distance Weight, IDW)空间插值方法进行空间分布特征分析。

2.4 趋势分析方法

Theil-Sen斜率(TSslope)估计是一种非参数估计法,用于估计时间序列数据的变化率,本文用于表征水体透明度的变化趋势。该估计具有处理删失回归模型的能力,并且对异常值不敏感。对于偏斜和异方差数据,可比非鲁棒简单线性回归更准确,即使对于正态分布的数据也能与非鲁棒最小二乘法竞争。尤其对于具有混沌特性的数据,其表达具有明显优势。Theil-Sen斜率的表达公式为：

(8)

式中,median为中位数函数,xj、xi为序列数据；tj、ti为序列数据对应的时间数据;序列长度为n、i、j为序号(1≤i≤j≤n)。当TSslope>0,表示上升趋势,反之为下降趋势；|TSslope|值越大则趋势越剧烈。

为了能够综合表达不同尺度湖泊水体透明度的变化率,本文提出了综合变化率,避免了单一维度的局限性。综合变化率为年均变化率(TS年)、季节(春、夏、秋、冬)变化率(TS季)及月(1—12月)变化率(TS月)的均值,如式(9)所示。

TS综合=mean(TS年+TS季+TS月)

(9)

其中mean为均值函数。

图3 技术路线图Fig.3 Model processing flow chart

图4 遥感影像MOD09GA波段之间的相关性分析 Fig.4 Correlation analysis among bands of MOD09GA remote sensing imagesGRA: 灰色关联度分析Grey Relational Analysis

3 结果与讨论

3.1 遥感影像数据的相关性分析

灰色关联度及Pearson相关系数分析结果如图4所示。由于波段之间的线性相关性较低,Pearson相关系数均处于较低水平,但可反映出波段间正负相关的性质,作为波段选择的辅助参考。由灰色关联度可得,b5/b6、1/b4、b3/b6、b2/b6、b1/b6、b2/b3和1/b5与水体透明度的相关性较高,能够较好的表达水体透明度的数据特征,所以选用此7个波段组合为模型构建的输入数据,对水体透明度进行反演。

3.2 水体透明度反演结果分析

图5 模型反演结果对比Fig.5 Comparison of model inversion results

以滇池2001—2010年的实测值与MOD09GA波段筛选后的数据为训练集与验证集构建水体透明度反演模型,并将此模型应用于2001—2018年滇池水体透明度的估算,模型评价参数如图5及表1所示,评价指标选用为均方根误差(Root Mean Squared Error, RMSE)和平均绝对误差(Mean Absolute Error, MAE)。对于具有混沌特性的数据而言,机器学习算法对数据特征的学习能力优于多项式模型,神经网络模型通过多隐含层多节点的学习,能够充分提取MOD09GA与水体透明度的特征。本文提出的估算方法(COD-LSTM)与其他模型对比结果表示,COD-LSTM模型估算效果最佳,且与实测值的趋势拟合较好(RMSE=0.1359, MAE=0.1134)。Binding等[34]仅考虑了单波段与水体透明度的对应关系,具有较大局限性；Wu等[9]和Knight等[19]同时考虑了两个波段与水体透明度的映射关系,拟合效果较优于Binding等[34]提出的方法,但对于水体透明度趋势的拟合具有较大局限；马建行等[18]提出的方法相对较优于以上三种方法,但水体透明度具有混沌特性,仅考虑多项式拟合较难挖掘数据特征,此方法仍较难准确估算水体透明度数值。本文基于深度神经网络提出的COD-LSTM模型,依靠遥感影像大数据的条件,以相关性分析、剔除异常值及去噪等方法充分过滤无效数据的干扰,并在多隐含层多节点的学习下,便于提取更多的数据特征,使得估算准确性得以大幅提高,因此在水体透明度的估算中具有更大优势。

表1 本文估算方法与其他反演方法及结果的对比

3.3 水体透明度的时空变化特征分析

图6 2001—2018年滇池水体透明度反演结果 Fig.6 The Secchi depth inversion results of Dianchi Lake from 2001 to 2018

2001—2018年滇池水体透明度的反演结果如图6所示,其中,综合变化率为年均变化率、季节(春、夏、秋、冬)变化率及月变化率的均值,这一指标能够综合表达不同尺度湖泊水体透明度的变化率,避免了单一维度的局限性。总体上看,滇池水体透明度呈下降趋势,综合变化率为-0.08 m/10 a,年均变化率为-0.05 m/10 a,季均变化率为-0.09 m/10 a,月均变化率为-0.08 m/10 a。按季节尺度分析,春季变化率最高,为-0.11 m/10 a,其次是秋季,为-0.09 m/10 a,夏季与冬季均为-0.08 m/10 a。按月分析,6月份变化率最高,为-0.15 m/10 a,最低变化率为-0.03 m/10 a,出现于10月份,高于月平均变化率的月份有5个(2—3月、5—6月、9月),低于月平均变化率的月份有6个(1月、4月、7—8月、10—11月),12月与月平均变化率持平。

2006年前后滇池水体透明度出现不同的变化趋势,2001年1月—2005年12月呈上升趋势,变化率为0.015 m/10 a；2006年1月—2018年12月呈下降趋势,变化率为-0.007 m/10 a,这可能与2003年开始实施的“一湖四片、南延北拓”的城镇化战略有关[35]。牛栏江-滇池补水工程于2008年12月30日开工建设,2013年底建成并投入试运行[36]。2009年1月—2013年12月呈下降趋势,变化率为-0.007 m/10 a；2014年1月—2018年12月呈上升趋势,变化率为0.001 m/10 a。由此表明,引水工程的实施对滇池水质的治理工作取得了一定成效。

2001—2018年滇池水体透明度年均值及综合变化率的空间分布如图7所示。水体透明度变化率的空间分布能够反映水体长时间尺度下透明度变化程度和趋势,能够为厘清透明度变化的空间异质性、有效保护治理水环境提供专题数据。研究结果表明,滇池水体透明度均呈下降趋势,水体透明度相对较高的区域,变化幅度相对较大,表明近18年来滇池水质总体呈不同程度的恶化趋势。滇池北部沿岸水体透明度较低,草海区域为全湖最低,但变化幅度较小,未表现出明显恶化趋势。南部和中部的水体透明度相对较高,但南部和中部的变化速率明显高于北部。滇池仅东北沿岸出现两个区域的水体透明度改善趋势,这可能与此区域建设的多个湿地公园有关。

水体透明度较低的区域分别位于主城区(滇池东部和北部)和海口镇(滇池西南部),变化速率较快的区域位于晋宁区(滇池南部),近年来,此区域的城镇化进程加快,可能是导致湖泊水质恶化的主要因素。滇池入水口位于北部,而北部恰好是城镇化扩张较为严重的区域,对流入滇池的水源的污染是巨大的。滇池水体透明度在2006—2013年间呈持续下降趋势,自2013年开始,昆明市政府投入约83.3亿元实施的牛栏江-滇池补水工程,但随着滇池流域“一湖四片、南延北拓”等城镇化战略的实施,滇池南岸和东岸近年来不透水表面扩张显著[23,27, 29],昆明市城市中心由1个变为了4个,2014年后水体透明度呈现明显的波动变化,表明高强度的城镇化过程导致的水质恶化抵消了牛栏江引水工程等水质保护和改善工程的效果,所以滇池水体透明度仍未有实质性的改观。

总体看来,2001—2018年间,滇池水体透明度呈下降趋势,从水体透明度时序变化特征能够反映出滇池水质仍处于恶化状态,尤其是距离城区较近的区域水体透明度相对较低,表明人类活动是造成水体透明度下降的主要因素；水体透明度较高的区域下降速度相对较快,可能是由于大部分区域的水质长期处于较差状态,在湖泊水动力过程的作用下导致状态较好的区域的水质发生了稳态转变。2014年后湖泊水体透明度下降速度得到了有效抑制,这与近年来政府部门的大力度监管与治理有直接关系。

图7 2001—2018年滇池水体平均透明度年均值及综合变化率的空间分布Fig.7 Spatial distribution of annual mean Secchi depth and comprehensive change rate of Dianchi Lake from 2001 to 2018

4 结论

水体透明度是评价水质标准的重要指标,但对较大空间尺度的监测具有一定困难,且由于研究区内针对湖泊水质原位监测工作起步较晚,导致湖泊水体透明度历史数据的缺失。为此,以MODIS遥感影像为数据,基于深度神经网络方法的长短期记忆神经网络模型(LSTM),对滇池2001年1月1日—2018年12月31日的水体透明度进行反演,并利用地理空间分析方法对其时空变化特征进行探究,填补了深度学习在水体透明度反演的空白,有效地提高了常用反演模型的精度,同时探讨了水体透明度的时空变化特征,为开展湖泊水质研究提供了新的思路和方法。

(1)本文基于深度学习方法提出的COD-LSTM模型在水体透明度反演中具有较好的性能(RMSE=0.1359, MAE=0.1134),能够充分提取遥感影像数据与水体透明度的特征,相对多项式拟合而言,具有较大优势。

(2)时序变化特征分析结果表明,滇池水体透明度总体呈下降趋势,综合变化率为-0.08 m/10 a,年均变化率为-0.05 m/10 a,季均变化率为-0.09 m/10 a,月均变化率为-0.08 m/10 a。2014—2018年间,滇池水体透明度以0.001 m/10 a的速率出现上升趋势,体现了牛栏江-滇池补水工程的治理成效。城市扩张的负反馈效应和牛栏江引水工程等环境治理工程的正反馈效应的叠加,抑制了水体透明度的下降速率,但水体透明度状况尚未发生质的改变。

(3)空间变化特征分析结果表明,水体透明度较高的区域下降率较大,水体透明度较低的区域变化趋势相对稳定,距离城区及居民区较近的水体透明度相对较低。人类活动将成为影响滇池水体透明度变化的重要因素,同时也是造成滇池水体污染的主要因素。随着中国城镇化的不断推进,人类活动对湖泊水体透明度下降的贡献将越来越大。

致谢：云南省环境科学研究院提供数据支撑,科技部国家遥感中心提供技术支持。