基于FUI模型的柬埔寨洞里萨湖水体颜色研究

许 杨，王 野，陆建忠*，冯 炼

(1.武汉大学测绘遥感信息工程国家重点实验室，武汉 430079；2.南方科技大学环境科学与工程学院，广东 深圳 518055)

伴随着人口增长而日益加剧的人类活动和全球变暖带来的气候变化，地球资源被过度消耗，生态环境也遭受到严峻的挑战，影响着全球生物的生存和发展[1-2].水环境是地球生态环境的重要组成部分，但近年来观察到很多水环境受到威胁的案例，如海洋溢油对沼泽的破坏、反复出现的大西洋马尾藻带等[3-4].从19世纪70年代海岸带水色扫描仪CZCS发射以来，人类可以利用卫星遥感数据更加宏观地监测海岸带水体[5-6].随着宽视场水色扫描仪SeaWiFS、MODIS的相继发射，因其在大范围监测水体环境的基础上，还能够更高精度地获取水体信息，可为全球湖泊、海岸带和海洋的水环境研究带来了新的机遇[7-11].SeaWiFS、MODIS等水色卫星数据在反演叶绿素浓度、总悬浮物浓度、透明度、溶解有机碳等参数方面的潜力也已经得到了有力的证明[12-15].2002年起中国相继发射了海洋水色系列卫星，其在水色、水温等海洋要素的反演及海岸带监测领域都有较好的应用前景[16-18].因此，在研究各类水体的水色要素和水质状态时存在着大量的卫星数据源可供选择，而具有更长观测历史的水色卫星数据可以挖掘出丰富的时空变化信息，能够为水环境的监测、治理和预报提供了更多的便利和应用价值.

已有大量的研究聚焦于使用卫星遥感数据来监测各类水环境要素，尤其以水色三要素(叶绿素浓度、总悬浮物浓度、黄色溶解有机物浓度)、固有光学量(吸收系数、散射系数、衰减系数等)、水华、富营养化程度等居多[8，19-25].近年来新提出了一种基于多光谱的传感器数据计算的水体颜色指数FUI，它是在国际照明委员会(CIE)色度系统定义的，源自海洋学中Forel-Ule (FU)观测.FUI用于刻画水体的颜色信息，并反映水质的综合状态[26].在色度空间定义的FUI不依赖于单一波段的信号响应，相比于水体光学成分的反演精度在其他组分浓度较高时会受限，水体颜色计算的不确定性可以被显著减小，使得FUI在水质状态的综合理解方面有着很好的应用前景[27].目前的研究已经证明了FUI与叶绿素浓度、水体富营养化程度、固有光学量等参数有较好的相关关系，大大降低了获取该类参数实测数据的成本[27-29].此外，海洋学领域中FU观测已经积累了上百年的历史数据，这为追溯水体光学参数的历史变化也带来了巨大的应用潜力[29].

洞里萨湖(Tonle Sap Lake)是东南亚最大的通河湖泊，是湄公河天然的蓄水池，对于调节当地生态系统和维持渔业资源等方面有着重要的意义.但是由于上游的建坝蓄水、滥砍滥伐和污染破坏，洞里萨湖的面积有减小趋势，水质状况堪忧[30].目前利用遥感数据针对洞里萨湖面积变化、湄公河三角洲的演变和洪峰强度等方面的研究已经有了大量进展[30-33].而对于洞里萨湖水色参数的研究多集中于总悬浮物浓度的反演及河湖交换的分析[32-35].因此，本研究将利用FUI模型来研究洞里萨湖水体颜色的时空分布，以期深入探究该湖泊水质综合状态及富营养化程度的变化规律，为洞里萨湖水环境的监测和治理提供重要的基础数据.

1 研究区概况与数据源

1.1 研究区概况

洞里萨湖(Tonle Sap Lake)位于柬埔寨境内的北部(11°53′～14°08′N，102°30′～105°29′E)，是东南亚最大的淡水湖泊.洞里萨湖在湄公河流域的下游，由洞里萨河与湄公河连通，其同湄公河的水流流向会出现季节性逆转的现象，属于季节性淹没的湖泊，具体位置及流域分布见图1.洞里萨湖属于热带季风气候，旱季在11月～5月，雨季在6月～10月[36].湖泊的面积和水量在一年之中会有明显变化，在旱季结束时，其面积和水量仅有2 500 km2和1 km3；而在9月和10月初，当湄公河水位达到最大值，西南季风降水也到达峰值时，洞里萨湖的面积和水量可分别达到16 000 km2和80 km3.近年来，由于全球变暖导致的气候变化及该区域人口增长带来的人类活动的影响，水体范围呈显著减小趋势，同时洞里萨湖的泥沙量增加，水质显著下降[30，37-39].

1.2 研究数据

本研究使用的数据是搭载在Terra卫星平台上的MODIS传感器三级数据——地表反射率产品MOD09A1(第六版)，数据来源于美国国家航空和航天局(NASA) (https://ladsweb.nascom.nasa.gov/).该数据是8 d合成的，空间分辨率是500 m×500 m，包括7个光谱波段，分别是蓝、绿、红、近红外和三个短波红外波段.在低级别数据通过大气校正，去除了气体、气溶胶和瑞利散射的影响后，还经过了严格的云和云阴影、气溶胶、观测角和蓝波段数值最小等筛选准则，最终每一幅影像都代表了该区域八天内的最佳观测影像[40].本研究共下载了866景2000年-2018年间覆盖洞里萨湖的影像数据.为了更加精确地刻画洞里萨湖的信息，在预处理阶段使用了MODIS重投影工具MRT，将所有的数据进行了重投影、裁剪.而后经过湖区范围内的去云筛选，共余下574景有效数据.

图1 研究区Fig.1 Study area

2 研究方法

本文拟采用王胜蕾等人提出的水体颜色指数FUI计算过程中的中间量α角，来研究洞里萨湖水色的时空变化[28].FUI指数是基于国际照明委员会CIE色度系统定义的，FUI取值1～21的离散值，能够实现从蓝色到棕色水体的全覆盖.而相应的色调角取36.64°～250.85°之间的连续值，能够更加精细地刻画水色及其变化，避免水体颜色及其变化的信息被压缩.

在使用MOD09A1产品计算色调角之前，需要应用一个简单的波段裁剪法来消除MOD09产品的噪声，并转为遥感反射率Rrs数据[41].接着使用红、绿、蓝波段的Rrs来计算CIE系统下的(X，Y，Z)坐标，转换公式为[42-43]：

X=2.768 9R+1.751 7G+1.130 2B，
Y=1.000 0R+4.590 7G+0.060 1B，
Z=0.000 0R+0.056 5G+5.593 4B.

(1)

再将新坐标归一化至0～1之间的值(x，y)，并使用公式(2)将其转换到图2所示的CIE色度图中[28].

(2)

此时，(x′，y′)处的向量和x′负轴之间的夹角定义为α角(以下简称色调角)，取值范围在0～360°.为了消除人眼感受到的真实色彩和传感器波段设置得到的颜色之间的误差，对色调角进行一个系统误差Δ校正.校正公式(3)如下：

Δ=-1.818 5×b5-62.164×b4+320.49×b3-

510.97×b2+306.32×b-50.343，R2=0.95,

(3)

其中，b=α/100.该公式是由辐射传输数值模型Hydrolight模拟的高光谱Rrs和MODIS等效波段之间建立的多项式模型，可以消除波段设置带来的误差[44].应用系统误差校正后，色调角可以重新定位到图2[28]中的位置，用于反映不同的颜色，而21个级别的水体颜色指数FUI是基于色调角而分段映射得到[43].考虑到FUI的取值范围(1～21)，进行系统误差Δ校正后的色调角有效值范围为36.64°～250.85°.

图2 FUI色调角及其在CIE色度图中的位置[28]Fig.2 The angle α and FUI shown in the CIE chromaticity diagram[28]

使用基于FUI的富营养化状态评估算法[28]，可以进一步研究多年来洞里萨湖富营养化程度的变化.该算法是根据Hydrolight仿真数据集与全球分布的实测数据集建立的，将水体分别归类于1) “贫营养”状态(1≤FUI<7)；2) “中营养”状态(7≤FUI<10，Rrs_645≤0.006 25)；3) “富营养”状态(10≤FUI<21，Rrs_645>0.006 25)[28].

此外，本研究还使用了归一化植被指数NDVI和归一化水体指数NDWI相结合的阈值法提取洞里萨湖的水面范围[20]，并将淹没频率大于80%的区域定义为常态湖区范围.使用该范围对所有色调角影像做掩膜.水面提取的具体步骤如下：1) 利用人机交互式图形用户界面GUI加载每景8 d反射率多波段影像，并且在GUI中生成相对应的NDVI (NDWI)和RGB真彩色图像，获取NDVI (NDWI)直方图，将两个峰的谷值作为初始阈值(Ti)；2) 通过初始阈值，获取水陆边界，将其叠加显示在RGB或NDVI(NDWI)图像上，如果边界与NDVI(NDWI)图像中的像素很好的吻合，则将Ti作为该副影像的最佳阈值.反之，手动调整Ti，直到获取最佳阈值和对应最佳水面提取范围[20].

如图3所示，以2001年2月18日影像为例，图3(a)～3(d)图依次是MOD09A1地表反射率产品合成的RGB影像、色调角影像、FUI影像和表征富营养化程度影像的示意图.该日的色调角取值范围在97.35°～228.64°，FUI取值范围为5～16，湖中兼有三种富营养化程度的水体.之后对每个月内的色调角、水体范围、FUI和Rrs_645数据分别计算一张中值图；再逐级统计该范围内季均、年均的色调角、面积和不同营养状态评级的面积比.

图3 研究区影像处理结果示意图(以2001年2月18日为例)Fig.3 Schematic diagram of image processing results in the study area (Take February 18，2001 as an example)

3 结果与分析

3.1 色调角的年际变化

本研究首先对色调角的年际变化进行了统计.如图4(a)所示，每个点即是色调角年均值，误差线表示该年色调角月均值的两倍标准差.19年间的色调角总均值是215.58°，FUI均值是14，即洞里萨湖的水体颜色一般处于从绿色到黄色的过渡色.通过统计发现，在2000年～2018年间，洞里萨湖色调角年际变化存在波动，整体上有略微下降的趋势，但不显著.其中，灰色方框中突出显示的2014年-2016年间的误差线，相比于前些年，明显变短.以上数据表明19年来洞里萨湖的色调角存在不显著的减小趋势，且2014年-2016年间各月水体颜色差异明显减小.Lin等人的研究表明洞里萨湖开阔水域面积的水文周期减小并呈现下降趋势，且季节性变化减弱.这些变化主要受到整个湄公河流域的水利大坝建设和相关农业灌溉的影响，而同期的气候影响关系不大[30].洞里萨湖面积的这一变化规律与本研究中色调角的年际变化趋势一致，因此本文将进一步探索色调角与面积的相关关系.

图4 2000年-2018年间色调角均值的变化趋势及与面积均值相关关系的统计分析Fig.4 Time trend of the average angle α from 2000 to 2018 and statistical analysis ofaverage angle α and average area

3.2 色调角的季节变化及与面积的关系

考虑到洞里萨流域雨季(6月～10月)、旱季(11月～5月)的典型划分，本研究分别统计了色调角月均值和面积月均值的时间趋势和相关关系，并进行了显著性检验，着重比较雨季、旱季不同月份的表现.图4(b)～(c)的结果表明，1) 2000年～2018年间洞里萨湖的色调角均值在不同的月份变化趋势不一致，旱季和雨季(蓝色背景)的月份均有呈现减小和增加趋势，但是由于云遮挡等原因，部分年份缺失有效月均值，呈现出的趋势并不完整.2) 面积均值在12个月均呈现为逐年下降的趋势，不同月份的色调角与面积表现出的相关关系正负不一致，旱季和雨季的月份均存在正相关和负相关关系.其中，处于旱季的12月达到最大的负相关关系(R=-0.86).处于雨季的6月达到最大的正相关关系(R=0.72).3) 色调角与面积相关关系的显著性检验结果显示，旱季和雨季均有4个月份的相关关系达到显著水平.2月、3月和8月处于显著水平(P<0.05)，而6月～7月、10月～12月的相关关系可被认为是非常显著的(P<0.01).所以无论洞里萨湖处于旱季还是雨季，多数月份的水体颜色变化很大程度上可以由面积变化来解释.

3.3 雨季、旱季典型月份的色调角变化率

在多数月份的色调角与面积呈现显著相关关系的背景下，研究选取了6月(雨季)、12月(旱季)这两个代表月份，进一步分析了2000年-2018年间洞里萨湖色调角或水体颜色的变化率及空间分布.图5(a)～(b)雨季代表月份的结果显示，湖中大部分区域的色调角在6月呈现出逐年减小或不变的趋势，少部分区域呈现增加的趋势.变化率被评定为显著和非常显著水平的区域集中于洞里萨湖北湖的湖中心偏东区域和南部的西岸，在6月这两个区域的色调角均呈现逐年减小的趋势.

图5 2000年-2018年色调角在雨季、旱季代表月份的变化率Fig.5 Change ratio of the hue angle in the representative months of rainy season and dry season from 2000 to 2018

图5(c)～(d)旱季代表月份的结果显示，湖中大部分区域的色调角在12月呈现逐年增加或不变的趋势，仅在洞里萨湖北湖的西岸沿线和南湖东岸的小部分区域，存在逐渐减小的趋势.面积相对更小的南湖基本都处于色调角逐年增加的趋势，且增加幅度高于北湖.12月份的显著性检验结果表明，逐年变化较强的区域多集中于北湖的岸边和南湖近岸区域.

总体上，在雨季的代表月份6月，洞里萨湖的色调角在大部分区域呈现出逐年减小的趋势；在旱季的代表月份12月，大部分区域的色调角则呈现相反的逐年增加趋势，且整体上逐年变化的幅度明显高于6月.无论在6月还是12月，色调角的变化率处于显著水平的面积较小，湖中大部分位置的逐年变化趋势并不显著.

3.4 雨季、旱季典型月份的各营养状态面积比

洞里萨湖存在着富营养化的问题，水质状态处于中等至较差水平，具有高营养程度和高固体颗粒物浓度的特点[45-47].本研究使用基于FUI的富营养化状态评估算法，获取了洞里萨湖2000年-2018年间不同营养状态的年均评级和各月评级.结果显示，各年年均评级均处于“富营养”状态，不同月份的营养状态评级存在差异，但没有统一特点.因此，本研究将着重关注6月和12月这两个雨季、旱季的代表月份，将基于19年间各营养状态所占面积比，分析洞里萨湖在这两个月份的营养状态评级及年度变化，得到如图6所示的结果.2000年至2018年间，洞里萨湖在雨季的代表月份6月，平均有99.97%的湖区面积被评级为“富营养”状态，其余被评级为“中营养”状态；在旱季的代表月份12月，平均有0.21%和4.65%的湖区被评级为“贫营养”和“中营养”状态，而被评级为“富营养”状态的湖泊面积约占95.14%.

在旱季12月，除了2000年、2001年、2011年和2013年外，其余各年都有90%以上的湖区处于“富营养”状态.在上述四年的12月份，分别有13%、20%、13%和15%的湖区处于较低的营养状态水平，即“中营养”和“贫营养”状态.现有研究中提及，2000年、2001年、2011年和2013年这四年的雨季洞里萨湖都发生了由于年降雨量异常导致的大洪水[48-49].此外，2002年该区域也发生了大洪水，但是该年洪水主导原因是上年的水位和水量较高[48].而本研究的结果显示，发生过由降雨量异常导致的大洪水后，在旱季12月洞里萨湖的营养状态似有明显好转的迹象，而其他年份的营养状态均以“富营养”状态为主.其他年份在旱季12月多以“富营养”状态为主，并且近5年没有明显的变化趋势，可见洞里萨湖在旱季的富营养化状态稳定且没有好转.

图6 2000年-2018年雨季、旱季代表月份不同营养状态面积比的变化Fig.6 Changes in the area ratio of differenttrophic statesin the representative months of rainy season and dry season from 2000 to 2018

4 结论与展望

本研究使用2000年-2018年间的MOD09A1地表反射率数据，在水体颜色指数FUI模型的中间量色调角的基础上，揭示了洞里萨湖在19年间水体颜色变化趋势和规律，并在同等统计规则下，分析了色调角与面积之间的关系，以及色调角变化率的空间分布和富营养化状态的评估.结果表明：1) 洞里萨湖在2000年-2018年间色调角呈现出下降的趋势，但不显著，湖泊的整体颜色由黄色向绿色转变.2014年-2016年间水体颜色的月间差异明显减小.2) 洞里萨湖各月的水体颜色变化趋势不同，19年间有7个月的色调角与面积均值为正相关关系，其中5个月份的相关关系呈现显著水平；有5个月份的色调角与面积均值为负相关关系，其中3个月份呈现显著水平，且这些月份在旱季、雨季均有分布.3) 就雨季、旱季的代表月份6月、12月的水体颜色变化率及空间分布而言，雨季6月的洞里萨湖以色调角逐年减小或不变的趋势为主，即水体颜色逐年偏绿；旱季12月则相反，水体颜色有变黄趋势，且变化幅度比雨季6月更加明显.在这两个月份中，湖中仅有少量位置处的水体颜色变化趋势达到显著水平，其空间分布也不统一.4) 19年间洞里萨湖在整体上呈现为“富营养”状态，在雨季的代表月份6月，平均99.97%的湖区处于富营养化状态.由年降雨量异常导致的大洪水年份中(2000年/2001年/2011年/2013年)，约有15%的湖区在旱季12月降为“中营养”或“贫营养”状态.其余年份的洞里萨湖在12月以“富营养”状态为主.

基于以上结论，本研究不仅能够为衡量洞里萨湖水体颜色及富营养化状态提供有效信息，并且可以在时间尺度上细化对其变化规律的研究.但是本文尚未涉及水体颜色变化规律的驱动因素分析，为了更加精确地量化人类活动和气候变化对洞里萨湖的影响，未来的研究可关注驱动因素与水体颜色参数、富营养化状态的关系.