基于光学和雷达多源遥感的于田绿洲土壤盐渍化时空分析

肖森天， 依力亚斯江·努尔麦麦提， 努尔比耶·穆合塔尔，赵 静， 阿迪莱·阿卜来提

（1.新疆大学地理与遥感科学学院，新疆 乌鲁木齐 830046；2.新疆大学绿洲生态重点实验室，新疆 乌鲁木齐 830046）

土壤盐渍化是全世界目前共同面临的生态环境问题，盐渍化问题会间接损害生态系统，严重阻碍干旱、半干旱区当地的农业生产活动和经济的可持续发展［1］。我国盐渍土壤分布广泛且各种盐碱土所占面积约9.9×107hm2，新疆的盐渍化土壤总面积约8.5×105hm2，是由当地气候条件和独特的地形地貌所造成［2-5］。利用遥感手段对绿洲地区土壤盐渍化时空分布和土壤盐渍化程度进行监测，实现当地自然生态环境稳定和确保绿洲农业生产的健康发展起到重要作用［6-8］。

遥感技术监测区域范围内的土壤盐渍化始于20世纪70年代［9］。国内外涉及遥感的应用中，光学图像的应用最为广泛，但由于光学遥感属于被动遥感，在获取数据过程中因其光谱特性易受到入射光源和天气条件影响，导致土壤信息的监测分类精度降低［10］，光学遥感利用地物的光谱特征对盐渍化土壤进行监测，但日间成像和受云雨影响的情况下，盐渍化信息的提取依赖光谱特征具有局限性［11］。

雷达遥感作为主动遥感，不受气候和昼夜影响，具备全天时、全天候的成像特点，且能提供光学遥感不包含的信息［12］。融合了微波遥感等多源遥感数据的监测方法可以弥补光学遥感数据的局限性，为目标检测提供光谱和纹理特征的独特组合［13］。目前基于微波遥感数据监测土壤盐渍化的相关研究中，学者们通过分析图像后向散射强度、土壤盐分和介电常数的关系，发现微波影像可以作为土壤盐分监测的有效手段［14］。利用雷达全极化数据，分析雷达后向散射系数与土壤盐分之间的关系，证明结合极化组合特征的微波影像能对不同程度的盐渍化进行区分［15］。使用光学和雷达影像多源数据，基于机器学习模型估算土壤盐分和空间分布，证明机器学习算法利用光学与雷达融合数据在土壤盐渍化监测方面的可行性［16-18］。

综上，集成多源光学与雷达数据应用在盐渍化监测的研究已有一定成果，但对盐渍化时空变化监测的相关研究较少。对此，本文探讨利用多源遥感数据和野外实测数据，通过提取有效的影像特征信息进行最优特征组合，对研究区使用随机森林分类方法，识别土壤盐渍化的分布并分析结果，研究于田绿洲土壤盐渍化的时空变化特征，探索微波与光学数据监测干旱区绿洲土壤盐碱化的能力。

1 研究区及数据来源

1.1 研究区概况

于田绿洲（36°30′～37°05′N，81°09′～82°03′E）位于塔克拉玛干沙漠南缘（图1），昆仑山北麓中部。东接民丰县，西邻策勒县，南与西藏自治区日土县相接，绿洲面积约3.95×104km2［19］。地势呈南高北低，垂直地带差异明显，相对高差4142 m，地形呈牛腿状［20］。绿洲内季节性河流发源于昆仑山的克里雅河属于以冰川融雪为补给的内陆河，流经于田绿洲，消失在塔克拉玛干沙漠腹地，当地水资源缺乏，水系相对封闭［21］。该区域属于暖温带内陆干旱荒漠气候，昼夜温差大，多年平均气温为12.2 ℃，最高温度39.4 ℃，最低温度-17.6 ℃；蒸发作用强，年降水量约33.5 mm；绿洲植被种类贫乏，生物群落结构单一，生态系统较脆弱，土壤盐渍化现象严重［22］。

1.2 数据资料与处理

1.2.1 遥感数据 雷达数据首先选取日本航天局发射的对地观测卫星（Advanced Land Observing Satel⁃lite，ALOS-2）所获取成像时间为2015 年4 月23 日的（Phased Array type L-band Synthetic Aperture Ra⁃dar，PALSAR）全极化数据。其次选择了欧洲空间局（European Space Agency，ESA）开发的Sentinel-1 卫星，成像日期为2018 年4 月23 日、2021 年4 月18 日的Sentinel-1A SAR 双极化数据。使用ENVI 5.3 软件中的SARscape 5.2.1 模块对影像进行多视处理、3 m×3 m窗口精细Lee滤波、地理编码、辐射定标，并重采样为20 m 分辨率，为与Sentinel-1 极化方式保持一致，选择VV、VH极化方式进行极化组合。

光学遥感数据采用Landsat 8陆地成像仪（Oper⁃ational Land Imager，OLI）数据，成像时间分别为2015年4月26日、2018年4月18日和2021年4月10日，具体参数信息见表1。使用ENVI 5.3 软件对时序Landsat 8 数据进行辐射定标、大气校正、影像配准、重采样为20 m分辨率等预处理流程。

表1 遥感数据参数Tab.1 Parameters of remote sensing data

1.2.2 野外数据 本研究根据于田绿洲的土壤盐渍化和交错带分布空间规律，于2015 年4 月22 日至5月7日和2021年7月9日至7月15日进行野外样品的采集，筛选出74 个有代表性的采样点，分布如图1a所示，每个采样点取0～20 cm表层土样、拍照并记录土壤类型、周围植被类型、建立采样点图片库便于对研究结果进行精度评价。将样本进行自然风干、研磨和使用1.5 mm 筛子过筛，使土壤样品和蒸馏水混合，溶液配制为1:5的土水比，使用电导仪在25 ℃的室温下测量电导率。通过构建电导率和总可溶性盐之间的方程式对土壤的含盐量进行计算和pH值的测定［23］。

图1 研究区影像及采样点Fig.1 Images and sampling points of study area

2 研究方法

首先获取于田绿洲研究区的野外实测数据、光学和雷达多源数据。然后将光学和雷达数据以多源融合的方式，通过特征提取、J-M距离分析和随机森林分类方法，对研究区土壤盐渍化进行分析，最后对研究区土壤盐渍化动态变化进行监测。具体技术路线如图2所示。

图2 技术路线图Fig.2 Workflow of study

2.1 特征提取及分析

遥感影像依赖不同的地物光谱信息区分地物类别，单一的特征信息无法提供足够的区分信息［24］。将不同的特征信息进行融合，融合后的特征数据在一定程度解决了单一特征信息片面问题，将各种特征描述的信息聚合，可提升高分辨率影像的分类精度［25］。

首先将预处理后的光学影像进行主成分分析（Principal Component Analysis，PCA），再对雷达影像进行后向散射特征提取并组合，最后使用Pearson相关分析，计算电导率EC值与后向散射组合特征间的显著性。根据图3 筛选出最优后向散射特征组合。以主成分分析后的光学影像为基础，将各分类特征进行叠加组合，建立分类特征数据集（表2），用同样的方法得到不同年份下的相同分类特征数据集。

表2 分类特征数据集组合Tab.2 Classification feature data set combination

图3 Pearson相关性矩阵Fig.3 Pearson correlation matrix diagram

2.2 J-M距离

为定量研究不同分类特征下地物间的可分离程度，对训练样本使用J-M 距离评估地物间的可分离性，J-M 距离的数值范围在0～2 之间，若两种地物的距离值越趋近2，说明两种类型的可分离性越高［26］。J-M计算方法如下：

式中：B表示巴氏距离。

计算结果如表3，方式三光学主成分特征和SAR 影像的极化组合特征叠加下，重度盐渍地与其他地物的J-M 距离最佳，表明各地物的可分离性在光学和雷达数据的组合下可区分性较好。

2.3 随机森林分类

将随机森林（RF）算法用于土壤盐渍化的分类。随机森林以集成分类器的形式进行分类，分类器由众多决策树组合而成，具有训练速度快、泛化性能强和模型稳定等优点［27］。随机森林方法对于集成多源遥感数据，具有多个分类特征的数据集可以有效提高分类精度［28］。徐乔等［29］对极化SAR 影像进行分类，利用随机森林分类算法实验得到高精度、光滑且分类层次连续的分类结果。本文结合文献［27-30］中的方法，根据特征数来对决策树的个数进行调试。

3 结果与分析

3.1 基于多源数据的盐渍化识别精度分析

为对比不同分类特征条件下地物的识别效果，将各年的PCA光谱特征影像和PCA+后向散射特征影像进行分类结果对比。因方式二，雷达影像进行多次分类后，总体精度平均为67.07%，难以达到区分地物的目的，故排除方式二。根据乔木等［31］对新疆土壤盐渍化程度分级标准，将地物分为裸地、植被、水体、轻中度盐渍地和重度盐渍地5 类（轻中度盐渍地为有盐生植被覆盖，pH值在7.1～9.5；重度盐渍地为无植被覆盖pH值在9.5以上）。本研究以高分辨率影像（奥维互动地图）和野外图片库作为目视解译的参考，选取训练样本和验证样本［32］。对不同的训练样本均使用相同的验证样本进行分类后的评价，便于科学对比分析分类后的验证结果。

根据表4 中分类精度可知，方式一的总体精度平均为80.36%，Kappa 系数平均为0.77。方式三的总体精度平均为85.62%，Kappa 系数平均为0.82。结合表3和表4，随机森林分类的精度说明多源数据特征组合的分类效果是要优于光学遥感数据。雷达影像的后向散射系数对土壤盐渍化的介电常数比较敏感，土壤含盐量直接影响介电常数的变化［33］，利用这一特性将雷达和光学影像结合，可以弥补雷达影像噪声对地物识别的干扰，增强后向散射系数对盐渍地的识别从而提高分类精度。

表3 典型地物的J-M距离Tab.3 J-M distance of typical ground objects

表4 分类精度验证Tab.4 Classification accuracy verification

为详细对比不同数据源分类的结果，以分类精度较好的2021年为例，选择于田绿洲具有代表性的Ⅰ、Ⅱ区域进行局部特征的对比。如图4所示，区域Ⅰ、Ⅱ是于田绿洲两个“牛角”边缘，沙漠与盐渍地交汇区域。这里土壤盐分较高，轻中度和重度盐渍地交错带十分明显。区域Ⅰ体现了绿洲、水体和盐渍地的过渡交错分布，PCA 光谱特征数据的分类结果中（图4a），两侧盐渍地存在错误分类的散点状植被区域，PCA+后向散射特征的分类结果较为有效的减少了零碎斑点的出现，体现盐渍化交错带的边界特征（图4b）。在区域Ⅱ地物由裸地沙漠、重度盐渍地和轻中度盐渍地逐渐过渡到水体区域，图4c和图4d 区域在植被和中轻度盐渍地的区分上存在略微的不同，因为植被比例在15%～70%，光学影像中地表反射率可能表现为植被和其他地物的混合［34］。根据野外采样发现耐盐植被如芦苇、柽柳和盐生草木可以在轻中度盐渍化土壤中生长，这些耐盐植被会间接影响光学影像中盐渍化土壤的光谱反射率［35］。PCA+后向散射特征影像具有的后向散射系数可以透过植被冠层，对地物的纹理信息进行捕捉，将其结合光学遥感的光谱信息，在一定程度上增强影像的分离度，消除植被的干扰，将轻度盐渍地的稀疏植被进行区分，降低其对地物分类的影响［36］。获取植被覆盖下的盐渍土含盐量信息，还需依靠微波散射模型来反演土壤含盐量［37］。因此，将被动的光学遥感和主动的微波遥感数据综合是目前监测土壤盐渍化较为理想的方法，未来研究中有必要结合土壤湿度、地表粗糙度和土壤介电常数来更好的解释土壤盐分的分布情况。

图4 随机森林分类结果对比Fig.4 Comparison chart of random forest classification results

3.2 土壤盐分空间动态变化分析

从分类结果可知（图5），于田绿洲盐渍化主要发生在北部、西南部，绿洲边缘和荒漠交错带区域，这与相关学者的研究一致［38-39］。盐渍化区域中，重度盐渍化土壤主要在研究区北部的两个“牛角”之间分布，以及克里雅河流域的下游区域，呈重度盐渍化向塔克拉玛干沙漠的过渡形式。在绿洲近郊轻中度盐渍化土壤分布在植被和重度盐渍化土壤之间，呈现条带状分布，主要原因是绿洲内排碱渠中的水未经处理排放到了绿洲边缘，导致绿洲边缘轻中度盐渍化程度较严重［40］。重度盐渍地与轻中度盐渍地分布区域交错分布，流经的区域因地下水位通常较高，伴随地表强烈的蒸发作用，土壤盐渍化积累作用显著［39］。

对研究区表层土壤盐渍化空间动态变化分析，提取各年同一时期盐渍化的面积变化（表5），于田绿洲2015 年的盐渍化面积占总面积的40.12%，其中，轻中度盐渍化区域占26.48%，重度盐渍化区域占13.65%。2018 年盐渍化面积占总体面积的35.18%，其中轻中度和重度盐渍化区域各占23.25%、11.93%。2021年盐渍化区域占总体面积的35.63%，其中轻中度和重度盐渍化区域各占20.51%、15.12%。2015—2021年于田绿洲的土壤盐渍化程度总体呈逐年波动下降趋势，盐渍化程度主要以轻中度盐渍地为主。轻中度盐渍地的面积变化从2015—2021年呈持续下降趋势，面积占比逐渐减少，植被面积呈上升趋势；重度盐渍地的面积呈先减少后增加的趋势；绿洲表层水体变化趋势不明显，但土壤水分是土壤盐分运动和变化的主要驱动因素，间接影响当地土壤盐分的空间格局变化［41］。

表5 2015—2021年地物分类面积变化Tab.5 The area change of classified features from 2015 to 2021

3.3 土壤盐分时间序列变化分析

为反映于田绿洲盐渍土在2015—2021 年的迁移变化，进行叠置分析得到盐渍地时空变化分布（图6），由于PALSAR-2影像范围的限制，2015年的研究仅在光学影像的分类基础上进行。使用桑基图直观可视化的体现了于田绿洲各地类之间的结构特征和转入转出关系（图7）。结果表明2015—2021 年于田绿洲的主要5 个乡镇区域的耕地周边，大部分区域处于盐渍化程度降低的变化趋势。但在喀尔克乡东南部的牧业村耕地存在盐渍化增加的现象，因绿洲内部土壤在排盐的改良措施下盐分含量较低，当灌溉和排放失衡时，土壤表层下的可溶性盐离子随毛管水运移到地表富积，造成耕地形成次生盐渍化［42-43］。盐渍化程度增加的希吾勒乡和斯也克乡北部区域，根据实地考察和DEM 影像，发现区域属于盆地地形，地下水的水位和矿化度相对较高。北部海拔高度相对于南部低约200 m，且干燥的气候导致水分蒸发，加重了土壤盐渍化的形成。

图6 2015—2021年盐渍地变化分布Fig.6 Distribution of saline change from 2015 to 2021

图7 2015—2021年地物类型转移图Fig.7 Feature type transfer map from 2015 to 2021

由研究区土壤盐分时间序列变化统计结果可知，2015—2021 年于田绿洲盐渍地面积变化量为-6723.31 hm2·a-1，变化率为-10.68%，说明盐渍化程度总体处于降低趋势。其中轻中度盐渍化年均变化量为-1526.58 hm2·a-1，变化率为-22.06%，呈降低的变化趋势；重度盐渍化的年均变化量为406.03 hm2·a-1，变化率为11.38%，重度盐渍化程度的变化过程是较为缓慢的增加趋势，说明重度盐渍化土壤的改良较为困难，这与前人对于田绿洲不同程度盐渍化土壤变化的预测结果基本一致［44］。

4 结论

以于田绿洲为研究对象，采用野外实测数据、光学数据和雷达数据结合的方式，通过特征提取和J-M距离分析，选择最佳分类特征组合，进行随机森林分类。得到以下主要结论：

（1）对比2015年、2018年和2021年分类结果和不同分类特征下采用随机森林方法对土壤盐分分类的精度，得到光学影像的PCA光谱特征的总体精度平均为80.36%，Kappa系数平均为0.77。PCA+后向散射特征的组合总体精度平均为85.62%，Kappa系数平均为0.82。精度评价结果表明使用随机森林分类方法在该分类特征下分类精度较高。说明在光学影像分类的基础上加入雷达影像是未来盐渍化监测中较为有效的方法手段。

（2）通过对于田绿洲土壤盐分空间动态变化分析，得出2015—2021年于田绿洲产生土壤盐渍化的区域主要分布于研究区北部的绿洲边缘和荒漠交错带。盐渍地的面积变化为轻中度盐渍地呈逐渐降低趋势，于田绿洲的植被面积逐年增加，重度盐渍地为先减少后增加的波动趋势，于田绿洲内的河流水体面积变化不显著。

（3）对于田绿洲土壤盐渍化时序变化分析，数据表明2015—2021年盐渍化程度总体呈下降趋势，盐渍地面积年均变化为-1120.55 hm2·a-1，变化率为-10.68%，说明近7 a土壤盐渍化现象得到较为有效控制和改良。