探究冰雪归一化消融指数的消融区自动提取.docx

探究冰雪归一化消融指数的消融区自动提取南极冰盖表面消融影响冰盖物质平衡和全球能量平衡[1]一方面,表面消融形成的融水以三种方式影响南极冰盖的物质平衡[2,3],融水汇集后从冰盖边缘流出[4];渗入裂隙后加剧裂隙的扩张,破坏冰架稳定性甚至导致冰架的大范围的快速崩解;渗入冰盖内部,促进冰盖内部通道的发育,抵达冰盖底部后减弱冰盖底部运动阻力,导致冰盖底部滑动[5,6]另一方面,南极地表97.8%为雪覆盖[7],干雪具有高反照率,而蓝冰、湿雪、融水反照率较低[8],其吸收的太阳辐射约为干雪吸收的三倍[9,10,11],形成“反照率低-吸收太阳辐射多-加剧融化”的正反馈机制[12,13],因此,消融范围及持续时间的变化将影响全球的辐射收支[5,14]尽管机制有所差异,但蓝冰、湿雪和融水都是冰雪消融的产物[3,15,16],本文统称它们为消融区在大多数年份南极冰盖表面消融区面积占比都小于20%,远低于占比常年超过40%的格陵兰冰盖,但南极冰盖表面融化的面积是格陵兰的两倍多[9]目前,格陵兰冰盖消融研究体系较完善[17,18],南极冰盖消融研究相对缺乏[1,19]南极冰盖消融研究更多为大尺度冰盖整体的研究,区域性精细研究较少且多集中于南极半岛[20]。

从消融的观测对象上看,冰面融池及河流水系的研究相对深入[21,22],而湿雪和蓝冰监测有待更多的开展遥感技术是冰盖表面冰雪消融大范围观测的最佳手段湿雪的范围及持续时间的变化可以通过微波辐射计、散射计、SAR以及光学数据监测前三者为微波数据不受云的影响而光学数据易受到云的影响微波辐射计数据对于消融的监测主要依赖于干湿雪亮温的显著差异[23,24],散射计数据对于消融的监测主要依赖于干、湿雪后向散射系数的显著差异[25],两种数据主要监测大尺度区域冰雪表面消融起止和持续时间、范围随时间的变化,无法满足对消融特征空间细节方面的描述[26,27,28]SAR数据对于湿雪的监测也主要依赖于干、湿雪后向散射系数的显著差异[29,30],其高空间分辨率可描述消融空间细节特征,可用于精确监测湿雪范围[31,32],但重访周期较长,适用于精确的年际变化分析,无法满足季节性分析[16,26]光学数据对于湿雪、融水及蓝冰的监测主要依赖于水、冰和雪光谱特征的差异例如利用融水和冰雪在蓝光和红光波段的光谱差异使用冰雪归一化水体指数开展冰面融水提取[4,22];基于光学数据对南北极水系水文特征的分析[4,22,33,34];基于Landsat-7影像数据蓝冰与湿雪裸岩的光谱反射率差异对整个南极蓝冰区的提取分析[35]。

光学数据用于湿雪的分析却鲜少开展综上所述,冰雪消融的监测难以同时满足高时、空间分辨率,高时间分辨率的被动微波数据难以满足消融细节特征的描述,高空间分辨率的SAR及光学数据难以实现季节性变化的相关分析,而两类数据很难融合因此,目前亟需高时、空分辨率冰盖消融的季节性分析本文选用Landsat-8数据开展东南极达尔克冰川消融区监测,其多光谱信息可用于湿雪、蓝冰和融水的提取,其较高的空间分辨率可满足对冰盖表面消融细节特征的监测此外,达尔克冰川刚好处于Landsat-8不同轨道影像的重叠区,重访周期为2—7天,可满足融化季节内多时间序列消融观测本文提出并验证了基于Landsat-8影像的MNDWIice消融区自动监测方法,以此为基础,对2016—2017年东南极达尔克冰川消融区的时空分布及季节性特征进行了初步分析1、研究区与数据1.1研究区概况达尔克冰川位于东南极拉斯曼丘陵东侧(69°26′S,76°27′E),邻接埃默里冰架(图1a),距离中山站直线距离不超过3km,为距离中山站最近的冰川达尔克冰川流域起源于南极冰盖内部(图1b),流入接地线东部的冰源为自上游600多千米的高原,而流入接地线西部的冰源为自上游仅不到10km的冰盖边缘(图1c),总面积超过1900km[2][36]。

达尔克冰川南部接地线距离冰川前缘约15km,漂浮部分即冰架面积约为20km[2](图1),是典型的冰盖溢出型冰川冰川前端入海口宽度为3km,随着冰川的不断运动而崩解[37]1.2研究数据本研究采用2016年9月1日—2017年4月30日覆盖达尔克冰川无云或少云的Landsat-8影像,此时间段内覆盖达尔克冰川的数据有41景,目视解译可用于研究的少云或无云影像共有18景影像(表1)数据来自USGS网站:https://glovis.usgs.gov影像由于采集时间的不同,太阳天顶角不同,日地距离也可能存在差异,为了降低这些差异的影响,使得不同时相影像数据具有可比性,将Landsat-8Level-1数据的像元灰度值(DN值)转换为大气表观反射率(TOA反射率)反射率计算采用Kang等[38]的方法,计算公式来自Landsat手册(http://Landsat.usgs.gov/),计算所需参数从每个Landsat-8数据集的头文件中读取此外,考虑到高纬度地区的同一幅影像,太阳高度角从中间到边缘的差异较大,计算采用逐像元计算的太阳高度角替代原始使用的单一中心点太阳高度角[38]下文中所做的计算及制图全部基于反射率数据。

2、研究方法2.1消融区的光谱特征利用多光谱光学数据识别目标地物的思路是利用目标地物与背景地物在不同波长上的波谱反射特性差异,选择最佳的波段进行比值运算方式,通过选择合适的阈值最大程度地区分目标地物与背景地物本文利用Landsat-8影像识别的消融区,包含裸露的蓝冰(无论是否发生融化)、湿雪和融水图1达尔克冰川概况.蓝冰生成过程中,冰内的气泡被挤出,使得冰晶变大,蓝冰呈现蓝色而湿雪和融水由于水的存在也显示蓝色,三者在遥感真彩色波段上均显示蓝色,区别于干雪(图2)Landsat-8影像不同地类的光谱信息如图3所示,蓝冰、湿雪和融水在近红外波段(波段5)表现出强烈的吸收而干雪表现出强烈的反射,在该波段消融区(蓝冰、湿雪和融水)与干雪、阴影和云差异最为明显,在蓝光波段(波段2)各地类差别最小基于此,本文选择利用近红外波段(波段5)和蓝光波段(波段2)进行比值运算,对消融区进行提取表1本研究阈值试验所用Landsat-8影像*其中有云指原始数据云量>20%,少云指0%<云量≤20%图2蓝冰、湿雪、融水、干雪、阴影样本对比图.底图皆为2016年12月21日Landsat-8反射率数据,显示方式皆为波段(5,3,2)对应(红绿蓝)标准差拉伸增强图像2.2改进的冰雪归一化消融指数本文采用的改进的冰雪归一化消融指数MNDWIice是基于归一化水体指数NDWI及冰雪归一化水体指数NDWIice提出的。

归一化水体指数NDWI是对绿光波段(G)与近红外波段(NIR)进行归一化差值处理,以凸显陆地表面的开阔水域[39],其表达式为:冰雪归一化水体指数(NDWIice)是对蓝光波段(B)与红光波段(R)进行归一化差值处理,以凸显冰盖表面融水[22],其表达式为:NDWIice最早提出是用于对格陵兰冰盖表面融水的提取[22],目前已用在南极表面融水深度的确定[4]本文基于消融区在Landsat-8数据的近红外波段(波段5)和蓝光波段(波段2)的光谱特征的强响应,提出改进的冰雪归一化消融指数MNDWIice为:式中,Bblue和Bnir分别为Landsat-8蓝光波段(波段2)和近红外波段(波段5)的反射率基于2016年12月21日获取的两期影像,选取融水、蓝冰、湿雪、干雪、云和阴影6类地物的纯像元样本,计算各地类对应NDWI、NDWIice和MNDWIice的指数值,如图3(b)所示试验显示MNDWIice对于消融区的区分效果最好,其他两种指数除对融水的区分度较明显外,其余几类地物的指数值相差较小,很难进行区分图3各地物类别反射率及指数特征.2.3基于F值法的阈值获取和消融区提取精度评估为选取最优阈值,引入机器学习评估方法——精度F值法(F-Measure),用F值评估分类结果精度。

F值由准确率(Precision,P)和召回率(Recall,R)两个度量值计算得到,F值越高表明提取效果越好准确率即为影像分类中的用户精度,召回率即为影像分类中的生产者精度(制图精度)目标像元分为四类:TP,TN,FP,FNTP是指消融区被正确划分为消融区的像元数,TN是指干雪被正确划分为干雪的像元数,FP是指干雪被错误划分为消融区的像元数,FN是指消融区被误判为干雪的像元数[40,41]准确率P用于评估被正确提取的消融像元占所有被提取为消融像元的比例,召回率R用于评估被正确提取的消融像元占应该被提取为消融像元的比例F值法是准确率和召回率的加权调和平均根据召回率与准确率对于研究的重要性给β赋值,当两者同样重要时,β值取1;当准确率更重要,β值取小于1的值;召回率更重要时,β值取大于1的值[40]本研究需同时考虑消融区被准确分类为消融区的概率和分类结果影像中消融区像元对应原始影像为消融像元的概率,因此定义β值取1根据MNDWIice阈值选取和分类结果精度评估需要,我们分别选取了试验样本和验证样本(如图1所示),为保证阈值选取结果的准确性,试验样本选择位于消融区和干雪区交界处的容易发生误分的区域,图1(d)中7个白框即为试验所用的试验样本区(包含干雪区及消融区),试验样本的面积占研究区面积的9.94%,各期影像的干雪及消融区试验样本即基于各期影像反射率数据对白框区域采用目视解译获取;验证样本选择的是30个空间随机分布的纯消融像元或者纯非消融区像元,每个验证样本像元大小为5×5,像元需满足纯样本条件,即像元区域为纯消融区像元或者纯非消融区像元,由于消融初期消融区面积小,20160907期影像难以选择一定数量的消融区纯样本,考虑到验证样本的一致性,故从20160914期影像选取30个纯消融样本和30个纯非消融样本,后续影像由于阴影等其他影响不满足纯像元条件的,选择其他纯像元做补充,保证验证样本为纯消融或纯非消融样本。

试验样本与验证样本存在一定比例的重合,但各期影像中重合面积不超过试验样本面积的0.95%,不超过验证样本面积的22.34%我们旨在利用单一的MNDWIice阈值进行不同时相的消融区的分类,单一阈值TMNDWIice计算方法如下:式中,TMNDWIice为18期影像利用F值法选取到的每一期最优阈值,Fmax(i)为最优阈值对应的F值每一期的最优阈值TMNDWIice获取方法:根据经验设定MNDWIice阈值范围为[0.01,0.2],以0.001为步长,利用试验样本数据计算每个阈值分类结果对应的F值,获取MNDWIice阈值-F值曲线,选择F值最大值Fmax(i)点处的阈值为最优阈值(图4)图420161221期影像阈值选取实验本文基于Landsat-8level1数据计算反射率,利用反射率数据计算得到各时相MNDWIice数据,利用18期影像的MNDWIice结合样本区目视解译结果确定各时相最优阈值及单一阈值,最后利用目视解译的验证样本结果,对消融区自动提取结果进行精度评估,具体实现流程详见图53、结果3.1基于MNDWIice消融区提取效果评估基于试验样本的阈值试验,18期影像消融区提取MNDWIice最优阈值如图6a所示,各时相MNDWIice阈值分布范围为0.1—0.17,根据每一期的最优阈值和对应的F值,依据公式7计算了18期影像消融区提取的单一阈值TMNDWIice为0.136。

采用各期影像最优阈值对其余18期影像试验样本的提取结果进行F值法精度评估(图6a),Fmax(i)值整体优于65.9%,准确率R优于56.4%,召回率P优于65.3%;采用单一阈值0.136对其余18期影像试验样本的提取结果进行F值法精度评估(图6b),Fmax(i)值整体优于63.。