梁林林, 江利明, 周志伟, 陈玉兴, 孙亚飞
(1.中国科学院测量与地球物理研究所大地测量与地球动力学国家重点实验室,武汉 430077; 2.中国科学院大学,北京 100049)
冻土的面积在北半球占比非常大,多年冻土和季节冻土面积分别占北半球陆地表面积的24%和30%[1-3]。中国是世界第三大冻土国,以青藏高原为主体的多年冻土面积约为1.5×106km2,约占我国多年冻土面积的70%[4]。全球气候变暖及过多的人类活动导致了地下冰融化、多年冻土退化等冻土问题,容易诱发热融滑塌,其形成的主要原因是在人为因素和气温升高的影响下,地下冰融化,融化土体产生不均匀沉降,覆盖的草皮、泥炭层及其下伏土层失去支撑而塌落下来[5]。热融滑塌一旦形成就会产生溯源侵蚀滑塌,发展速度快、滑塌范围大,对冻土区工程施工以及工程建设后的长期稳定性影响非常大[6],同时对周围地区环境也会造成难以自然恢复的破坏,如加速了多年冻土退化,导致冻土地区水土流失、植被退化,加剧冻土荒漠化等。大范围获取冻土区域热融滑塌分布状况及其时空变化信息,对冻土区工程建设规划以及冻土灾害防治具有指导意义,同时也可用于冻土退化研究和冻土水碳释放定量估算,对研究多年冻土退化与全球气候变化相互影响具有重要意义。
对高空间分辨率影像使用传统的基于像元的分类方法会产生“椒盐”现象,分类精度低,不能很好地利用影像的空间信息。而面向对象分析(object-oriented analysis,OOA)是以影像对象作为分类的最小单元,通过对对象进行分析,训练提取对象的光谱、形状和纹理等多种语义信息,使用分类方法实现遥感影像的高精度分类,弥补了传统的基于像元分类方法存在的语义不足,更适合高空间分辨率影像。随着遥感和计算机技术的发展,OOA技术已在遥感影像地理空间对象提取中被广泛应用。Daniel等[7]通过对SPOT5和数字高程模型数据应用基于面向对象的半自动分类方法对意大利西北部的滑坡进行提取,通过预先存在的滑坡库存数据和PSInSAR结果进行验证,证明了该方法的有效性; 甘甜等[8]利用玉树县结古镇震后的QuickBird高空间分辨率遥感影像,采用OOA技术实现研究区内建筑物震害信息提取,为灾害损失快速评估提供了可靠的基础数据; 孙娜等[9]利用高分一号PMS影像,采用OOA的影像分析技术与区域种子生长方法相结合实现了黄土高原地区水体的高精度提取; 王岩等[10]利用汶川地震震后高空间分辨率航空遥感数据对建筑物震害情况进行了面向对象的快速提取与自动分类,震害提取结果显著,表明OOA技术在地震研究中具有可行性。
卫星遥感影像因其空间和时间分辨率较低而难以对研究对象进行精细监测,并且大气和云层的干扰也可能使得数据无法使用。近年来,无人机在环境监测中的应用得到了快速发展,它可以按照要求获取高空间和高时间分辨率的影像,弥补了遥感卫星影像的缺陷。无人机技术在遥感监测调查中具有无可比拟的优势,在自然科学领域中得到了越来越多的使用[11-13]。
本文利用Trimble UX5无人机获取了青藏高原黑河流域上游东支的俄博岭垭口冻土区高空间分辨率遥感影像,研究了面向对象的冻土热融滑塌边界检测方法,详细分析了遥感影像分割和分类等关键技术,利用最邻近、K-最邻近(K-nearest neighbor,KNN)、决策树、支持向量机(support vector machine,SVM)和随机森林5种面向对象分类方法对冻土区热融滑塌进行识别检测,为研究冻土时空变化和冻土水碳储量估算提供了关键基础资料。
研究区位于祁连山黑河上游东支的俄博岭垭口(E100°55′,N38°00′)区域(图1),海拔为3 685 m,处于甘肃省和青海省的交界处,丝绸之路G227国道旁,著名历史古城峨堡镇附近。研究区内年平均气温低于4 ℃,年平均降水量为350 mm,大陆性气候明显。研究区多年冻土发育,因而具有典型的湿地生态系统特征,植被类型为高寒草甸。
图1 研究区影像
由于全球气候变暖的影响,冻土年平均地温升高,地下冰和积雪融化,雨水下渗,研究区内多年冻土退化,发育着典型的热融滑塌,如图2所示。研究区内存在2条典型的冻土热融滑塌区域,大滑塌区域长度约400 m,面积约7 600 m2; 小滑塌区域长度约120 m,面积约1 200 m2。
(a) 热融滑塌区域1 (b) 热融滑塌区域2
采用天宝公司生产的UX5固定翼无人机(图3),搭载Sony A5100数码相机作为遥感影像采集平台。影像获取时间为2016年7月28日,航拍飞行高度为200 m,航向重叠度和旁向重叠度都为80%,共获取高空间分辨率影像96幅,包括红光(R)、绿光(G)、蓝光(B)3波段,空间分辨率为5.22 cm,单幅影像覆盖面积为0.06 km2。本文所使用的其他数据包括由R,G,B波段组成的数字正射影像(digital orthophoto map, DOM)和数字表面模型(digital surface model, DSM),以及由DSM影像计算得到的坡度(slope)数据。其中DOM和DSM均由UASmaster软件处理无人机航拍单幅影像生成,主要包括空中三角测量、点云、DSM和DOM影像生成等步骤。DOM影像的空间分辨率为0.06 m,DSM和坡度数据的空间分辨率均为0.14 m。数据覆盖范围大约为0.34 km2。
图3 UX5固定翼无人机
野外实测数据是利用天宝R8型RTK对研究区内2条典型的热融滑塌边界进行野外采集所得,图4为全球定位系统(global positioning system,GPS)现场测量。野外采集使用的是连续地形测量的方法,每隔0.5 m左右采集一个点数据,每个点测量的时间为1 s,围绕滑塌边界测量一周,共测量1 800个点,经处理形成2条热融滑塌的边界,用于热融滑塌提取精度评定。
图4 GPS现场测量
本文主要采用eCognition软件对无人机DOM,DSM和由DSM数据得到的坡度数据通过最邻近、KNN、决策树、SVM和随机森林等5种不同的面向对象分类方法提取热融滑塌信息。首先,对坡度、DOM和DSM数据进行多尺度分割,并分别选择不同的分类特征组合; 然后,采用5种面向对象的监督机器学习方法完成热融滑塌信息提取; 最后,对提取结果进行精度评价与分析。
在进行面向对象分类时首先要对影像进行分割,将影像的同质区域分割出来,作为面向对象分类的最小单元,本文使用多尺度分割算法来实现影像分割[14-15]。该算法是一种自下而上的分割算法,从任一个像元开始,根据相关的异质性标准与其相邻的对象进行区域合并,通过迭代形成小于设定同质性尺度阈值(即分割尺度)的影像对象[16]。异质性标准由影像对象的光谱异质性和形状异质性相结合计算而得[17],其中形状异质性包括光滑度和紧致度2种属性。异质性标准f计算公式为
f=ω1hcolor+(1-ω1)hshape,
(1)
式中:ω1为光谱信息权重,取值在[0,1]之间;hcolor为光谱异质性;hshape为形状异质性。
光谱异质性hcolor为特定权重值下各个数据层光谱值的标准差之和,计算公式为
(2)
式中:pi为第i波段的权重;σi为第i波段对象的标准差;n为波段数。
形状异质性hshape计算公式为
hshape=ω2u+(1-ω2)v,
(3)
式中:ω2为紧致度权重;u为影像区域整体紧致度;v为影像区域边界光滑度。
影像区域紧致度u及光滑度v计算公式分别为
(4)
v=E/L,
(5)
式中:E为影像区域实际的边界长度;N为影像区域的像元总数;L为影像区域最小外接矩形的边界长度。
影像分割过程中,需要考虑上文提到的波段权重、光谱信息权重、紧致度权重和分割尺度等参数。波段权重是衡量各个波段在数据处理中的重要性,因为本文使用的数据各波段影响相同,故设置各波段权重相同,权重值为1; 光谱信息权重是衡量分割过程中光谱信息相对于形状所占的权重; 紧致度用于根据较小的差异区分紧凑目标和不紧凑的目标; 分割尺度则决定分割对象的大小,既能够将地物边界显示清楚,又不能太破碎,本文为了比较分割参数对分类结果的影响,通过反复实验,确定了分割参数,如表1所示。
表1 分割参数
面向对象分类以分割对象作为分类的最小单元,分析提取对象的纹理、光谱和形状等信息,并将这些信息作为先验知识加入到分类器,通过一定的分类方法实现高空间分辨率遥感影像的高精度分类。为了验证该方法的可行性,本文利用最邻近法、KNN、决策树、SVM和随机森林等5种面向对象监督机器学习方法对研究区内的热融滑塌进行提取。
监督机器学习方法首先需要选择训练样本进行训练学习,样本选择的好坏直接影响分类结果。训练样本的选择遵循了均匀分布和完整性2大原则。尽量在研究区内广泛而均匀地为待分类类别分别选择样本,使选择的样本具有代表性和典型性,既能够很好地反映地物的光谱分布,又能减少“同物异谱”和“异物同谱”现象。
光谱特征是遥感影像中最直观的信息,能够反映出不同地物的特性。但由于无人机高空间分辨率影像的波段数较少,光谱信息比较有限,只利用光谱特征会存在一些不足。由于研究区内热融滑塌区域的坡度值变化会比较明显,因此选择了地形信息来提取热融滑塌。本文对研究区内的光谱信息、几何信息和地形信息进行了综合研究,分割尺度和分类特征如表2所示。为了实验分割尺度和分类特征对实验结果的影响,分别基于相同分类特征不同分割尺度和相同分割尺度不同分类特征采用上文提到的5种面向对象监督机器学习方法来检测研究区内热融滑塌(分别记为SY1,SY2和SY3),实现对研究区内热融滑塌分布检测。
表2 分类特征提取参数
本文利用最邻近法、KNN、决策树、SVM和随机森林等5种面向对象监督机器学习方法,按照表2的分割参数和分类特征对研究区进行了3次热融滑塌边界检测实验。以SY1实验为例进行结果展示,其中图5(a)为样本分布,图5(b)—(f)为各方法提取结果及GPS实测数据,为清晰展示,仅显示热融滑塌类的分类结果。
(a) SY1样本分布 (b) 最邻近分类结果 (c) KNN分类结果
(d) 决策树分类结果 (e) SVM分类结果 (f) 随机森林分类结果
为了客观评价不同分类方法的分类结果,验证面向对象分类方法应用于热融滑塌边界检测的可行性与准确性,本文以实地测量的GPS数据作为真实滑塌边界,按照滑塌与非滑塌区域面积比例在ArcGIS软件中随机选取83 063个样本点作为验证数据建立混淆矩阵对分类实验进行精度评价,分类精度评价指标如表3—5所示。
表3 SY1分类混淆矩阵
表4 SY2分类混淆矩阵
表5 SY3分类混淆矩阵
分析表3—5可知,3次实验总体精度均高于90%,分类效果较好;制图精度均高于73%,最高达到94%,表明热融滑塌边界识别效果良好,有利于在没有先验知识的情况下发现冻土热融滑塌问题。SY3的用户精度明显低于SY1和SY2,主要因为分类特征较少,不能很好地区分滑塌与非滑塌区域;而SY1和SY2的用户精度相差较小。
5种分类方法的Kappa系数如图6所示。分析图6可知,SY1和SY2的Kappa系数基本在0.6以上,最高达到0.72,且2次实验Kappa系数相差较小,而SY3的Kappa系数则明显低于SY1和SY2,最低只有0.35。表明在本次实验中,分割尺度对分类精度影响较小,分类特征对分类精度影响较大,因此选取合适的分类特征是提高分类精度的关键。5种方法中,3次实验SVM的Kappa系数均高于其他4种分类方法,分别是0.72,0.69和0.61,表明SVM方法比较适合该研究区热融滑塌边界的检测提取。
图6 5种分类方法Kappa系数比较
本文利用无人机高空间分辨率影像和实地测量GPS数据,采用5种面向对象分类方法对俄博岭垭口冻土区热融滑塌边界进行了检测及验证,实验结果表明:
1)采用无人机高空间分辨率影像和实地测量GPS数据,利用面向对象方法对研究区热融滑塌进行提取,提取总体精度均在90%以上。表明基于无人机高空间分辨率影像的面向对象方法对提取青藏高原冻土热融滑塌边界的可行性和准确性,有良好的应用发展前景。
2)对比3次分类实验结果,分割尺度不同而分类特征相同(SY1和SY2)的Kappa系数相差甚微,表明本实验分割尺度对分类精度影响较小; 分割尺度相同而分类特征不同(SY2和SY3)的Kappa系数相差比较大,SY3的分类效果较差,表明分类特征对分类精度影响较大,因此分类需要选择合适的分类特征。
3)对比分析3组分类实验的Kappa系数,SVM方法的Kappa系数高于其他4种分类方法,表明基于SVM面向对象监督学习的无人机影像分类方法更适合于该研究区热融滑塌边界提取。
志谢:感谢香港中文大学刘琳教授和兰州大学张廷军教授的指导。