基于Sentinel-1多时相InSAR影像的云南松切梢小蠹危害程度监测

薛娟， 俞琳锋， 林起楠， 刘广， 黄华国

(1.北京林业大学林学院省部共建森林培育与保护教育部重点实验室，北京 100083；2.中国科学院遥感与数字地球研究所，北京 100094)

0 引言

森林病虫害是威胁森林健康的重要因素。近年来，国内外森林病虫害呈增长趋势，严重程度也在不断增加，合理有效开展病虫害监测对森林保护具有重要意义。森林病虫害发生面积大，传统的实地调查耗费人力物力，监测难度大。随着科学技术的发展，遥感作为新兴的监测手段被广泛应用，可以在大尺度范围快速高效地获取森林病虫害信息[1]。光学遥感依据植被的光谱特征进行监测[2]，特别是随着高光谱技术的发展更能实现对森林病虫害的有效监测[3-4]。但光学遥感存在“异物同谱”问题，且易受云雨天气影响而造成数据的缺失，在多云雨地区的应用受到限制[5]。合成孔径雷达(synthetic aperture Radar，SAR)具有全天时、全天候获取数据的特点，同时可提取后向散射、干涉和极化等信息[6]，可为森林病虫害信息提取提供更有利的数据支持。

云南松是云南地区的主要树种，近年来受切梢小蠹危害严重。切梢小蠹是既蛀干又蛀梢的害虫，有全年危害、生活隐蔽、繁殖快速的习性，易对当地树种造成区域性危害，给林业发展带来巨大损失。森林病虫害的大面积爆发以及云南地区多云雨的气候条件限制了传统的监测方法和光学遥感的应用。SAR数据的后向散射系数与干涉相干系数受多种因素影响[7]，其中林分的生长状态、环境因素、受害类型、土壤和枝叶密度的差异都会导致后向散射系数和相干系数的差异。基于这些差异特性，本文提出利用SAR数据开展森林病虫害的监测方法研究，并以云南省祥云县为例，基于对多时相Sentinel-1数据的后向散射系数、干涉相干系数与云南松物候及相对湿度的相关性分析，同时通过对多时相数据的融合，成功实现了采用多时相SAR干涉数据对森林病虫害的分类。该方法可以有效区分多云多雨条件下健康林、轻度受害林和重度受害林。

1 研究区概况

研究区位于云南省大理州祥云县(N25°12′～25°52′，E100°25′～101°02′)。该区地形以山地为主，亚热带高原季风气候，年降雨分布不均，主要集中在6—9月，雨季降雨量和云量较大，9月—次年5月间则干旱少雨，年均降雨量为810.8 mm，相对湿度为68%，年均气温为14.7℃。这样的气候条件使当地易形成大面积的持续干旱。祥云县森林覆盖率高达65%，其中以云南松纯林为主，森林稳定性较低。近年来切梢小蠹的大面积爆发导致云南松林遭受了严重破坏，给当地林业发展带来巨大损失。

图1 研究区地理位置

切梢小蠹属鞘翅目小蠹科，是一种严重危害松属树种的蛀梢且蛀干的害虫，在云南地区一年繁殖一代，极具隐蔽性，多数时间隐藏在树体内部[8]，全年均取食寄主树木，危害树木生长。其危害包括蛀干期和蛀梢期2个阶段：成虫在树干上侵蛀繁殖，卵期、幼虫期和蛹期在树干内危害；成虫羽化转到枝梢上发育，汲取枝梢营养物质[9]。

本研究试验样地选择在普淜镇天峰山云南松林区。通过对样地内单木的总梢与健康梢进行计数得到单木枯梢率，所有枯梢率的算术平均值即为该林分的枯梢率，进而利用枯梢率对云南松林进行切梢小蠹的危害程度划分。根据林业有害生物发生(危害)程度标准，并结合当地的受害情况，将枯梢率低于10%的林分划为健康林，10%～50%的林分为轻度受害林，50%以上的林分为重度受害林[10]。

2 研究数据与方法

2.1 遥感数据

本研究选用欧空局的哨兵-1(Sentinel-1)数据。该卫星是欧洲航天局哥白尼计划中的地球观测卫星，于2014年4月发射升空，传感器为C波段雷达，可以获取包括条带成像(stripmap,SM)、干涉宽幅(interferometric wide swath,IW)、超宽幅(extra wide swath,EW)和波浪(wave,WV)4种成像模式的数据，具备单极化(HH/VV)或双极化(HH+HV/VH+VV)的数据获取能力，单颗卫星重复周期为12 d，最高空间分辨率为1.73 m×4.3 m(距离向×方位向)[11-12]。本研究使用Sentinel-1A的IW模式VV极化数据，获取时间从2015年5月—2016年9月，包括19景升轨数据和12景降轨数据。此外，采用SRTM DEM作为辅助数据支持。

2.2 其他数据

2.2.1 野外采样点

2016年5月对研究区进行了大面积野外调查，大致确定了松林受害程度，并在每类危害程度的松林中选择代表性样地，目测数出枯梢个数和总数，计算得到枯梢率。根据林地的生长情况、枯梢率以及受病虫害严重性均匀布置30个30 m×30 m的样地，进行每木检尺。其中健康林样地8个、轻度受害林样地10个、重度受害林样地12个；同时对研究区进行GPS定点，估测样本点30 m以内的枯梢率；并结合Landsat8 OLI影像进行目视解译，增加样本点共55个，用于分类结果的验证。

2.2.2 气象数据

气象数据来源于Reliable Prognosis(http: //rp5.ru)，由气象数据国际自由交换系统提供，选择祥城镇站点的数据。

2.2.3 物候数据

对云南松的针叶生长期、树高生长期、开花结实期和球果生长期等定性指标进行2015—2016年2 a间的连续观测，并根据观测指标进行物候期的划分。

2.3 数据预处理

运用SARscape对Sentinel-1数据进行预处理，包括影像配准、辐射定标和多视化处理、相干性估计、滤波和地理编码等。由于SAR数据的成像特点，地形起伏会对影像造成很大影响，需要去除地形影响[6]，所以基于90 m空间分辨率的SRTM DEM数据进行地形校正，最终得到空间分辨率为20 m的时间序列的强度影像和干涉影像。对影像取3像素×3像素窗口的均值作为样地的影像数值。具体处理过程如下：

1)强度影像处理。对多时相数据的影像进行配准，对于升轨、降轨数据分别以2015年5月16日和5月30日为基准影像在SARscape中进行自动配准，达到亚像元配准精度；选择多时相雷达滤波器进行滤波处理，过滤斑点噪声；基于DEM影像对SAR数据进行地形地理编码，同时进行辐射定标和叠掩/阴影处理，获得后向散射系数(图2(a))。

2)干涉影像处理。估算影像间的基线情况，结果表明所有影像的时间基线和空间基线均在临界范围内，可以进行干涉处理；利用DEM数据进行配准，去除平地效应，生成去平后的干涉图；对去平后的影像进行自适应滤波处理，进行相干性估计，生成干涉图；利用DEM数据赋予影像地理坐标系，消除地形畸变，生成地理编码的相干系数(图2(b))。

(a) 后向散射系数 (b) 相干系数

2.4 相关性分析

地物的后向散射系数与相干系数受多种因素影响，包括地物表面粗糙度、体散射系数、介电常数等环境目标参数和雷达波波长、极化方式、入射角度等传感器参数[13]。本文主要通过相关分析的方法来研究病虫害危害程度与后向散射系数、相干系数之间的关系。

1)单景相关性分析。用单因素方差分析法分析单景图像后向散射系数和相干系数与病虫害危害程度相关性的显著性[14]。

2)时序相关性分析。单个时相的健康林与受害林差异可能不够明显，因此将全年影像全部加入，绘制后向散射系数和相干系数随时间变化的趋势图，进行危害程度及趋势的分类评价，并融合多时相Sentinel-1数据进行健康林与受害林的分类。

3)进行物候期与气象数据因素贡献分析，用于解释时序变化。

4)结合实地野外调查数据与Landsat8 OLI影像目视解译数据进行验证。

3 结果和讨论

3.1 单景相关分析

单景影像方差分析及多重比较结果如图3所示(显著性水平α=0.05)。

(a) 升轨数据相干系数 (b) 降轨数据相干系数

(c) 升轨数据后向散射系数 (d) 降轨数据后向散射系数

从图3可以看出，相干系数(图3(a)和(b))中，有12景升轨影像与10景降轨影像的p值小于0.05，但在多重比较中，这些影像均存在p值大于0.05，说明在α= 0.05水平下，单景影像的健康林、轻度与重度受害林总体存在差异，但多重比较中，差异性不显著；后向散射系数中(图3(c)和(d))，有9景升轨影像与3景降轨影像的p值小于0.05，但多重比较中，这些影像的轻度与重度受害林的p值均大于0.05，说明在α= 0.05水平下，健康林与受害林总体存在差异，但轻度与重度受害林不存在差异。

综上所述，单景影像的相干系数与后向散射系数均不能完全区分健康林、轻度受害林以及重度受害林。其中单景的后向散射影像不能区分轻度与重度受害林；相干系数中有单景影像可以区分轻度与重度受害林。可能的原因是重度受害林中存在大量死木和倒木，是重要的散射体，所以后向散射系数中受害林没有很大的差别，而轻度与重度受害林中健康枝梢的比例不一样，所以相干系数有差别。

3.2 时序相关性分析

时序分析结果如图4所示。

(a) 升轨数据相干系数 (b) 降轨数据相干系数

(c) 升轨数据后向散射系数 (d) 降轨数据后向散射系数

从图4可以看出，相干系数Coh(图4(a)和(b))中，健康林、轻度与重度受害林的时序变化趋势一致，且具有明显的季节性，夏季的相干系数低于冬季，这与夏季气候因素变化大有关；与升轨数据相比，降轨数据的相干系数中重度受害林、轻度受害林、健康林的相干系数依次增大，健康林与受害林的区分度较大，但轻度与重度受害林的差别较小；雷达干涉的相干系数与植被生长的变化、森林冠层结构以及地表情况有关，也会受到风、降雨等气候因素的干扰，在健康林中地表情况比较稳定，因此相干系数较高，而受害林中存在大量灌木，这些灌木林具有较强的去相干作用，使受害林的相干值降低。

后向散射系数γ°(图4(c)和(d))中，健康林、轻度与重度受害林的季节性变化均不明显；升轨数据中重度受害林、轻度受害林、健康林的后向散射系数依次增大，但健康林与受害林的区分度较小；对于针叶林而言，雷达波段可以在一定程度上穿透其冠层，其后向散射主要来源于树干的反射和地表的二面角散射，但由于雷达采用斜视成像的模式，其后向散射系数在山区受地形的影响大，所以健康林与受害林的后向散射系数差别较小。

时序相关性分析表明利用相干系数进行病虫害程度分类，健康林与受害林的可分离度较高，特别是降轨数据。但是，轻度与重度受害林的相干系数区分度较小；利用后向散射系数进行分类，健康林与受害林的差别较小，分离度较低。

3.3 物候分析

根据云南松针叶生长期、树高生长期、开花结实期和球果生长期等定性观测指标，将云南松物候期划分为抽梢期(4—5月)、针叶生长期(6—10月)、休眠期(11月—次年2月)。具体物候期生长特点见表1。基于云南松物候期对相干系数与后向散射系数的时序变化进行分析。图4(a)与(b)分析表明，相干系数与云南松物候期相关，整体趋势是“低—高—低”。云南松在休眠期较稳定，相干系数整体高于生长期。从云南松物候期看，4—5月和7月份是转折点，相干系数也是最低的；4月份开始进入抽梢期，云南松新梢生长，针叶也开始抽长变绿，云南松林变化较大，所以相干系数变低；7月份进入针叶生长盛期，是相干系数较低的重要原因之一；之后高生长和粗生长变缓，云南松逐渐进入休眠期，相干系数增高。而图4(c)与(d)分析表明，多时相影像间的差别没有相干系数大，即后向散射系数与云南松物候期相关性较低。综上所述，相干系数能够指示云南松物候期各个生长节点。

表1 云南松物候期及生长特点

3.4 气象数据分析

相对湿度与林分含水量有一定联系，可能与微波散射信号有关。基于地面高度2 m处的相对湿度分析相干系数与后向散射系数的时序变化。通过分别对相干系数(Coh)、后向散射系数(γ°)与相对湿度进行相关性分析(*表示在0.05水平上显著相关，**表示在0.01水平上显著相关)，如图5所示

(A为升轨数据，D为降轨数据)。

图5 相对湿度与SAR参数的相关性分析

对图5分析，相干系数与相对湿度基本不相关。相比之下，后向散射系数与相对湿度有一定的相关性，其中受害林的相关性比健康林高，轻度受害林的后向散射系数与相对湿度相关系数最高，升轨数据在0.05水平显著相关，降轨数据在0.01水平显著相关，相关系数达到0.78。Sentinel-1 SAR数据为C波段，受害林与健康林相比，树冠较稀疏，所以健康林能更好地反映树冠表层信息，而受害林林分的环境因子相关性较高。

3.5 验证结果

将获取的样本分为训练样本和验证样本，其中分类样本30个、验证样本55个，使用ENVI软件进行监督分类处理。分类结果见图6，分类精度结果见表2。

(a) 升轨数据相干系数 (b) 降轨数据相干系数 (c) 升轨数据后向散射系数(d) 降轨数据后向散射系数

表2 分类精度验证

从分类结果和精度评价看，健康林与受害林成片分布，这与云南地区切梢小蠹成片危害的实际情况相符。利用相干系数进行分类精度较高，特别是降轨数据分类总精度达到83.15%，其中健康林和重度受害林分类精度均较高；轻度受害林分类精度较低，易被错分为健康林或重度受害林，枯梢率接近健康林与重度受害林的林分分类精度低。后向散射系数分类总精度相对较低，特别是降轨数据受害林的分类精度都很低。

4 结论

云南多云雨的气候条件极大限制了光学数据的应用。本文利用多时相C波段InSAR数据基于时间序列的干涉影像和后向散射强度影像以及云南松物候、相对湿度资料，分析后向散射系数与相干系数的时序特征，开展森林病虫害程度分类，并基于地面数据对融合多时相数据进行了分类精度验证。主要结论如下：

1)后向散射系数和相干系数的时序变化与云南松物候相关，且相干系数的相关性较高，采用相干系数能够指示云南松物候期各个生长节点。

2)相干系数与相对湿度基本不相关。相比之下，后向散射系数与相对湿度有一定的相关性，其中受害林的相关性比健康林较高，特别是降轨数据的相关性达到0.78。健康林能更好反映树冠表层信息，而受害林与林分的环境因子相关性较高。

3)单时相SAR数据中，健康林与受害林的后向散射系数与相干系数差别很小，需要综合SAR数据的时相信息与干涉信息，将多时相InSAR数据进行融合参与分类。实地精度验证表明多时相InSAR可用于森林病虫害监测。