基于NDVI滇池水华特征的时空变化研究

何云玲，熊巧利，罗贤 ，李同艳，余岚

1. 云南大学资源环境与地球科学学院，云南 昆明 650091；2. 云南大学亚洲国际河流中心，云南 昆明 650091；3. 云南省国际河流与跨境生态安全重点实验室，云南 昆明 650091

水质是水体的物理（色度、悬浮物等）、化学（有机和无机物含量）和生物（细菌、微生物、浮游生物等）方面的特征及其组成状况。生活和工农业生产废水中大量氮、磷等营养物质进入淡水水体后导致藻类（主要由蓝藻、绿藻组成）大量繁殖的“水华”现象，是湖泊富营养化的典型特征。水华一直是评价湖泊水体水质状况的重要指标之一。自上个世纪 90年代以来，中国太湖、巢湖、滇池等湖泊蓝藻水华暴发相继出现，对饮用水源安全造成较大影响（李旭文等，2011）。在采取一系列治理措施后，“三湖”水质有所改善，富营养化程度较历史峰值均有一定降低，然而以蓝藻为优势类群的水华现象却依然时有发生（Zhang et al.，2015；Zhang et al.，2016）。如何监测和控制藻类水华，是湖泊富营养化研究中的重点和难点（杨桂山等，2010；蒋大林，2015）。

水质遥感是指利用遥感技术对水质指标进行检测，并根据遥感波段的信息与水质指标的光谱特征，建立反演模型对水质指标的浓度进行预测（冼翠玲等，2017）。相比于传统的水质监测手段，遥感监测技术具有扫描面积大、资料更新周期短、成本相对较低和便于进行长期实时监测等优点。自20世纪 70年代初期开始，国外就广泛地开展水质遥感反演研究（Song et al.，2012）。国内的相关研究在最近几年发展迅速，闻建光等（2006）基于Hyperion星载高光谱遥感影像，运用波段比值法、一阶微分处理技术构建了太湖叶绿素a浓度的反演模型。Shi et al.（2015）发现MODIS Aqua影像与总悬浮物在645 nm波长处具有极强相关性，据此发展了适应性极强的总悬浮物经验反演模型，并在太湖取得了良好效果。归一化植被指数（Normalized Difference Vegetation Index，NDVI），是遥感影像中近红外波段的反射值与红光波段的反射值之差与两者之和的比值，常被用于确定目标区绿色植物覆盖程度的指标值（陈旭等，2015）；徐京萍等（2008）发现利用MODIS数据的NDVI能够确定湖泊蓝藻的时空分布；周立国等（2008）基于MODIS卫星影像，采用比值模型方法确定蓝藻分布的相对浓度信息；蓝藻暴发时引发水体温度、色度和透明度等一系列物理性质发生变化，进而导致水体反射波谱特性也发生变化，利用这一特性，就可以基于NDVI值对蓝藻水华分布范围和面积进行识别（鲁韦坤等，2009）。

滇池位于云南省会昆明市的南部，是该市重要的水源地，支持着昆明市工业、农业、旅游业等发展，历史上素有“高原明珠”之称。但随着社会经济的快速发展，大量的工农业废水和生活污水排入滇池，加之入湖支流水系较多，出流水系较少，湖泊换水周期较长，致使滇池水质恶化，水体富营养化问题突出（王杰等，2018），成为全国湖泊富营养化治理的重点对象。已有学者对滇池水质的时空分布特征做过研究，主要集中在滇池水质传统方法采样监测（梁中耀等，2014）、叶绿素 a浓度分析（韦力元等，2016）、NDVI值域的划分（谢国清等，2010）、水质等级预测模型（徐玉妃等，2018）、以及气象条件与水质之间的关系方面（张虎才等，2017）。本研究基于2000－2017年的NDVI不同值域对滇池蓝藻水华信息进行提取，并从逐年、逐月到逐旬详细分析滇池水华的时空分布规律以及变化特征，以期为滇池的治理及水资源的开发利用提供数据支持，也为建立实时监测和预警系统奠定基础。

1 滇池水域概况

滇池，亦称昆明湖，是云南省最大的淡水湖，位于昆明市主城区的西南方向，湖面海拔1886 m，平均水深5 m，最深8 m；湖水在西南海口泄出，称螳螂川，为长江上游干流金沙江支流普渡河上源，属长江流域金沙江水系，为地震断层陷落型的湖泊。南北长39 km，东西最宽为13 km，流域面积约为2920 km2，容水量为15.7亿立方米。滇池流域地处低纬度、高海拔地区，属于热带高原季风气候；年均温 14.7 ℃，日温差较大，降雨年内分配特点明显，分为干湿分明的雨季和旱季。流域内水系众多，有 30余条河流呈向心状注入滇池，其中径流面积大于100 km2的河流有盘龙江、宝象河、洛龙河、捞鱼河、大河、柴河、东大河。滇池流域的地理位置及基本概况如图1所示。

2 数据来源与处理

2.1 数据来源

图1 滇池流域地理位置及基本概况图Fig. 1 Geographical location and basic profile of Dianchi Basin

利用中国科学院资源环境科学数据中心（http：//www.resdc.cn/）2015年中国土地利用遥感监测数据提取滇池水体边界，分辨率为30 m×30 m。DEM 数据来自于美国国防部国家测绘局公布的SRTM-DEM 数据（http：//www.cgiar-csi.org），分辨率为90 m×90 m。MODIS数据来自于美国国家航空航天局 NASA（National Aeronautics and Space Administration）（https：//ladsweb.modaps.eosdis.nasa.gov/)的Terra遥感卫星2000年2月18日－2017年12月31日的MOD13Q1数据产品，选用16 d合成250 m×250 m分辨率的NDVI数据产品；该数据产品已经完成影像数据的辐射定标、大气校正、几何校正、气溶胶、臭氧吸收订正及去云等预处理。

2.2 数据处理

将MODIS影像数据导入MRT软件中进行影像拼接并进行UTM投影转换。借助ENVI 5.1、ArcGIS 10.2软件进行空间裁剪、掩膜、重分类得到 2000－2017年的时间序列遥感数据集。并采用最大值合成法（maximum value composite，MVC）建立NDVI时间序列栅格数据集，分别求出月、季、年 NDVI值，公式如下：式中，i为月序号，取值范围1－12月；NDVIi为第 i月 NDVI最大值；NDVIk为月上下旬 NDVI值，k取值为1－24；j为季序号，取值范围为1－4，NDVIj为第j季NDVI最大值；y为年序号，取值范围为1－18，NDVIy为第y年NDVI最大值，yNDVIj为第y年第j季NDVI最大值。

采用最小二乘法计算数据集中所有像元的NDVI与时间的回归斜率Slope，若Slope＞0，表示随时间变化NDVI呈增加趋势，并且数值越大增幅越大；若Slope＜0则代表NDVI呈下降趋势。计算公式为：

式中，i为年序号；n为年跨度，即18年；NDVIi为第i年NDVI值。

采用谢国清等（2010）在滇池蓝藻水华光谱特征研究中关于NDVI阈值划分水华的等级方法，即NDVI值大于-0.1的区域，为蓝藻水华区域，其中-0.1≤NDVI≤0.2为轻度水华，表征像元内水华覆盖度为0－30%；0.2＜NDVI≤0.4为中度水华，表征像元内蓝藻水华覆盖度为 31%－80%；NDVI＞0.4为重度水华，表征像元内蓝藻水华覆盖度为 81%－100%。应用GetData Graph Digitizer 2.26软件对以往实测研究中滇池富营养化指数变化特征进行数值提取，对比分析及验证基于遥感识别基础上NDVI表征的滇池水华特征。

3 研究结果与分析

3.1 年际变化

从2000－2017年滇池NDVI的逐年变化趋势图中可以看出（图2），通过实测研究得到的2008－2014年滇池富营养化指数与同时段基于 NDVI表征的水华特征变化趋势相比，出现较高值和较低值的波动变化趋势相似。其次，除了秋季有缓慢增加的趋势外，其他季节和全年的NDVI整体上均呈减少趋势，但均未通过P＜0.05的显著性检验。据统计，冬季的NDVI减少趋势最明显（-0.0028/a），夏季次之（-0.0022/a），春季变化趋势不明显。另外，夏季NDVI与年NDVI的变化趋势较为一致，18年间NDVI经历几次增加和减少的波动。总体而言，从2000年开始滇池年NDVI呈减少趋势，蓝藻水华覆盖度呈降低趋势。

图2 2000－2017年滇池NDVI值与富营养化指数变化图Fig. 2 Annual NDVI and eutrophication index in Dichi from 2000 to 2017

通过数据分析得到 2000－2017年滇池 NDVI的年平均值空间等级分布及变化斜率（图 3），多年平均的轻度水华面积占滇池总面积的54.99%，中度和重度水华面积分别占总面积的 27.42%和17.59%，从南到北，滇池水华呈现从轻度到中度、重度的变化格局（图3a）。滇池大部分区域NDVI变化斜率主要分布在 0－-0.1之间，NDVI整体上呈现不同程度的减小趋势（图3b），呈减小趋势的区域占研究区总面积55.64%；但有区域差异，NDVI减小，表征蓝藻水华覆盖度减少的地区主要分布在滇池北岸的断桥，滇池中部的观音山西、白渔口、海口西等地；NDVI增大，表征蓝藻水华覆盖度增加的地区主要分布在滇池北岸的草海、罗家营，滇池东南部的观音山东、滇池南等地。

图3 2000－2017年滇池年平均NDVI空间分布（a）及变化斜率（b）图Fig. 3 Spatial distribution of average NDVI and variation slope in Dianchi from 2000 to 2017

3.2 四季变化

按照上述 NDVI阈值划分蓝藻水华等级的方法，统计各季节滇池水华的空间分布特征，仍然表现出“北重南轻”的分布态势（图 4）；春、夏、秋三季滇池以轻度水华为主，分别占滇池总面积的86.49%、72.11%和 83.23%；冬季则以无水华区域所占总面积比例最大（69.34%），冬季的重度水华区域仅占1.57%；而中度和重度水华则主要出现在夏季，所占比例最高，分别占滇池总面积的13.89%和14%。总体而言，滇池的重度水华区域主要分布在滇池的北部（断桥、草海）和滇池的边界；中度水华区域主要分布在滇池的中部罗家营和观音山西；轻度水华区域主要集中在滇池南部的观音山中、海口西和滇池南。

图4 不同季节滇池水华特征的空间变化Fig. 4 Spatial distribution of water bloom in four seasons in Dianchi

3.3 月变化

通过实测得到的滇池逐月富营养化指数与基于NDVI表征的逐月水华特征变化趋势相比较（图5），可以看出二者出现较大值与较小值的波动变化相吻合，均是在7月达到最大值，12月达到最小值；5－8月滇池富营养化指数为高值区；同时，NDVI也表现为高值区，代表蓝藻水华覆盖度大，富营养化程度严重。

图5 滇池逐月NDVI值与富营养化指数变化图Fig. 5 Monthly NDVI value and eutrophication index in Dianchi

再从逐月的滇池NDVI值空间变化特征可以看出（图6），每年的1－3月除了靠近滇池北岸附近区域以外的大部分地区NDVI值很小，从4月开始，NDVI开始逐渐增大，在不同地区出现水华；7－8月份NDVI值较大，滇池水华较严重；9月份开始滇池NDVI值又逐渐减小，水华面积缩减。

统计不同程度水华出现的面积比例在各月的变化（图7），滇池11－4月以无水华区域为主，面积比例占绝大优势（尤其11－3月占80%以上）；5－10月则以轻度水华为主（尤其5－9月所占面积比例为80%左右）（图7a）。中度水华主要出现在5－10月之间，所占总面积比例均小于12%；重度水华则出现在 7－9月之间，所占总面积比例较小（均在8%以下），尤其重度水华出现的时间从5月开始到8月达到峰值（月平均值为7.61%），10月开始逐渐消失（图7b）。

3.4 旬变化

从上述滇池水华特征的四季和月变化分析中可知，蓝藻水华的暴发主要集中在夏季，进一步对夏季各旬的NDVI值空间变化进行研究（图8）。蓝藻水华的暴发时间7月上旬和8月下旬，其重度水华面积分别高达11.4%和14.0%。水华较严重的区域仍集中在滇池的草海、断桥和滇池的边界地区。

4 讨论与结论

4.1 讨论

滇池在维持流域生态平衡、满足生产生活用水、减轻洪涝灾害和提供丰富水产品等方面发挥着不可替代的作用。20世纪60年代其水质为Ⅱ类，水产资源丰富，70年代为Ⅲ类，80年代为Ⅴ类，90年代为劣Ⅴ类，“九五”开始，滇池被列为全国重点治理湖泊。伴随遥感监测技术的进步和迅速发展，对滇池的相关研究也逐渐开展，谢国清等（2010）利用假彩色合成法以及 NDVI阈值进行了滇池蓝藻水华划分，但未通过时间序列NDVI值对水华时空变化特征进行研究；王杰等（2018）通过传统采样方法对滇池 2008－2014年富营养化程度进行分析，但是数据样本较少；蒋大林（2015）通过浮游藻指数（FAI）对滇池水华时空分布特征及其驱动因子进行研究，主要集中在FAI阈值提取和水华面积的简单统计方面。本研究基于遥感数据，借鉴上述学者研究成果，对NDVI值域进行了分类，对18年来滇池水华的时空分布特征进行了更深入、详细的分析研究，研究结果和侯秀丽等（2018）、余丽燕等（2016）、余佑金等（2017）的结论相一致，表明基于NDVI值域判断滇池水华时空变化特征具有一定可行性。

滇池地处昆明城市下游，是昆明盆地最低凹地带，接纳生活污水、工业废水和含有农药化肥的农业污水，造成研究区蓝藻水华暴发格局呈现“北重南轻”的分布态势。此外，昆明地区 7－8月气温较高，风速较小，气候条件加速了蓝藻水华的暴发，导致了滇池水华在月变化上存在显著差异。2000年开始，昆明市开始大力推进滇池污染治理的各项工作，通过铺设截污管道、建污水处理厂、牛栏江-滇池补水工程通水，成效逐步显现，水质逐渐提升。这和本研究中2000－2017年滇池年平均NDVI总体呈现减少趋势，蓝藻水华覆盖度有所降低的结论是一致的。滇池北部的草海、断桥、晖弯和观音山地区过去是滇池流域富营养化程度较重、水质较差的地区，但是经过整治，这些地区的水质有变好趋势，反而是过去富营养化程度较轻、水质较好的滇池南部罗家营南端部分地区水华有加重趋势，这在未来滇池水质综合防治中值得关注。

图6 滇池水华的月变化空间分布图Fig.6 Spatial distribution of monthly water bloom in Dianchi

虽然基于MODIS遥感数据研究水华变化有其局限性，水体识别主要依赖于光谱信息导致NDVI无法区分蓝藻区域和高浑浊水体，不能有效识别低密度蓝藻水体，研究结果易扩大或缩小蓝藻的分布范围（郭鹏等，2012）。然而，不可否认的是，MODIS数据以其波段窄、信噪比高、波谱敏感性好以及可免费获取等优点可为实时、准确、经济地获取蓝藻水华信息提供便捷。在今后的研究中需探索更精确混合像元估算蓝藻水华面积的方法；需充分利用不同遥感数据的优点，提高数据时空分辨率。不同学者对NDVI值域划分有不同的划分依据，如谢国清等（2010）认为NVDI≥-0.1为蓝藻水华区域，孙佳丽等（2010）认为NDVI＞0为蓝藻水华区域。本研究采用的是文献引用量达高频次的分级标准，在今后的研究中可进一步分析通过 NDVI值域变化判断藻类水华的时空分布特征；以及进一步分析水华暴发与当地气象条件、经济社会相关因子的关系（王菁晗等，2018；鲁韦坤等，2017；欧阳志宏等，2015；盛海燕等，2015）。目前基于MODIS数据建立相应的蓝藻预警预测模型已有相关报道（魏清宇等，2008），随着卫星遥感技术的发展，传感器性能的进一步提高，直接通过遥感手段开展湖泊水质机理层面的分析以及预测预警，将是今后重要的研究方向。

图7 滇池轻度和无水华区域（a）中度和重度水华区域（b）所占面积比例图Fig. 7 Monthly variation of water bloom area ratio in Dianchi

图8 滇池水华的旬变化空间分布图Fig. 8 Spatial distribution of water bloom in every half month in Dianchi

4.2 结论

本研究基于滇池 2000－2017年时间序列的遥感影像数据，利用 GIS空间分析技术，对其近 18年来 NDVI的时间变化和空间分布特征进行了分析，探讨了蓝藻水华不同等级的分布面积和空间变化特征，得出以下结论：

（1）2000－2017年滇池NDVI总体呈现减少趋势（趋势系数未达到P＜0.05的显著性检验），呈减小趋势的区域占研究区总面积55.64%，蓝藻水华覆盖度有所降低。

（2）7－8月滇池NDVI值最大，水华最严重，9月开始水华面积缩减，11－4月以无水华区域为主（面积比例 80%左右），5－10月以轻度水华为主（80%左右），重度水华出现在7－9月（8%以下）。蓝藻水华的暴发时间为7月上旬和8月下旬，其重度水华面积分别高达11.4%和14.0%。

（3）滇池蓝藻水华的暴发格局总体呈现“北重南轻”的分布态势，水华较严重的区域集中在北部草海、断桥、晖弯和滇池边界。

（4）滇池NDVI减小表征蓝藻水华覆盖度减少，富营养化程度减轻的地区主要分布在滇池北岸的断桥，中部的观音山西、白渔口、海口西等地；NDVI增大表征蓝藻水华覆盖度增加，主要分布在滇池北岸的草海、罗家营，滇池东南部的观音山东、观音山中、滇池南等地。