基于高分一号卫星影像的延安市顾屯流域耕地盐渍化定量反演

李志刚，许 强 ，赵宽耀，陈婉琳，王晓晨，方汕澳，李为乐

(成都理工大学地质灾害防治与地质环境保护国家重点实验室，成都 610059)

黄土高原地区由于其特殊地质与气候环境，流域内水土流失较为严重，从而会引起耕地资源不足等问题。位于陕西省延安市宝塔区的顾屯流域就是一个典型的由于黄土的侵蚀造成的耕地减少的例子[1-3]。为解决土地资源需求问题，政府于2013年开始进行治沟造地工程，经过5年的治理，治沟造地工程已经基本完成，其中仅延安地区经治沟造地工程可增加的耕地面积就达330km2[4]。治沟工程均位于流域沟谷内，因而在形成新耕地时也改变了沟谷、土壤结构及局部的水文条件，环境改变可能会引发次生灾害问题，如滑坡、盐渍化等[5]。顾屯流域新增耕地由于水文条件改变、降雨量较少、蒸腾作用强烈等原因，部分区域有盐渍化加剧的趋势[3]。盐分的积聚会使土地退化，已造耕地难以使用，折减治沟造地工程的效益[6-8]。因此，有必要对顾屯流域的沟内耕地进行盐渍化定量评估，查明盐渍化的发展状况，为治沟造地工程的建设与耕地资源的保护提供指导。

目前已有学者对顾屯的盐渍化发展状况开展了一定程度研究。邹锡云等[3]对沟内库坝典型位置的盐分进行了提取，揭示了水库与盐渍化的形成有关。Jin等[4]对顾屯流域的主沟水库附近盐分信息进行了收集，发现流域上游相较于下游盐渍化程度较高。现阶段对于沟内的盐渍化研究都集中于对点状信息的提取，而点状信息难以反映整个沟内的盐渍化发展现状。遥感手段相较于实地测量，具有覆盖广、成本低等优点，对于区域尺度盐渍化的信息评估，往往采用遥感手段[9-11]。中外学者已经开展了很多利用遥感影像反演盐渍化信息的相关研究。江红南等[12]在对土壤信息和遥感信息分析的基础上提取了多个表征指数作为变量，建立了盐渍化信息提取模型，揭示了新疆渭干河-库车河三角洲土壤盐分特征。丁建丽等[13]选取土壤、植被和光谱特征为研究对象，筛选出对盐渍化最敏感的数据形式，对盐分信息进行了提取。冯雪力等[14]收集Radarsat-2雷达影像，通过多元回归建模对内蒙古河套灌区盐分分布进行了高精度监测。Zhang等[15]通过野外调查获取土壤盐分数据，分析Landsat 8影像与地表光谱特征，建立土壤盐度信息提取模型与植被指数模型。Wang等[16]利用Landsat 8的影像测试了不同波长差值、比值的线性回归模型，对艾比湖湿地盐渍化进行了估算。已有研究主要利用回归分析等方法对经验模型进行拟合，缺乏对有效波段和波段组合的筛选，并且传统的回归分析也难以解决有高精度需求的非线性问题。此外，顾屯耕地区面积相对较小，使用中等分辨率多光谱影像难以满足精确反演的要求。

现以顾屯流域作为研究区，利用2m高空间分辨率的高分一号(Gaofen-1)影像，首先利用多个传统光谱指数对耕地盐分含量进行定量反演；进而通过实测土壤样本对高分一号4个波段的信息量、波段相关性、波段与盐分的相关性进行计算，筛选特征波段组合，通过支持向量机(support vector machine，SVM)的机器学习方法挖掘特征波段与盐分信息之间的非线性关系，建立高精度的顾屯流域耕地盐渍化反演模型。该方法有助于治沟造地工程后耕地土壤盐分监测，对耕地资源管理保护有一定指导意义。

1 研究区概况

研究区为陕西省延安市的治沟造地工程区之一的顾屯流域(图 1)，位于宝塔区的甘谷驿镇，属半干旱大陆性季风气候，春夏时节较为干燥，雨季一般是每年的7—9月。当地年平均气温为10.3℃，年平均降雨量527 mm(2010—2017年)，年潜在蒸散发约1005mm，几乎是年降水量的两倍[4]。顾屯流域主沟总面积约为24km2，河道长度约12.5km，主沟坡度较为平缓，呈U形断面，两侧支沟较为陡立。沟道内土壤成分主要为黑垆土、黄绵土、红土等土壤，其中黄绵土占主要成分，具有较高的孔隙度与渗透率[17]。流域从2012年进行治沟造地后，通过开挖流域的边坡来填埋扩大沟底，从而形成耕地，同时修建水库和库坝，阻止土壤流失，显著改变了当地的物质及水分运移特征。在工程后形成耕地中，发现表面局部有盐分出露[图 1(b)、图1(c)]。

2 数据与方法

2.1 研究数据及预处理

2.1.1 遥感影像与实地采样

收集了研究区2018年6月11日的高分一号影像，影像全色波段分辨率为2m，另有分辨率为8m的蓝、绿、红及近红外4个多光谱波段。基于ENVI5.3软件对该高分影像进行预处理，包括辐射定标、大气校正及正射校正。对应影像时间段，在研究区盐渍化程度不同、空间分布均匀的56个位置进行现场采样[图 1(a)]，各采样点均匀取表层20cm深的土壤。

图 1 研究区位置图及对应位置天地图影像Fig.1 Location map of study area and images of the corresponding position

2.1.2 土壤样本处理与盐渍化分级

采用较常用的水土比5︰1浸提法测定收集土壤样本的盐分含量，具体步骤为称取过1mm筛的风干样本20g，加入蒸馏水100mL制得清壳滤液，进而测得盐分。对所有土壤的全盐量进行测定，得到有效样本56个。将样本按盐量由大到小排列，以4为间隔选取14个样本作为后续反演的验证样本，剩余42个样本作为建模样本。

目前盐渍化的等级划分没有固定标准。研究表明，土壤中全盐量在1g/kg时对农作物的生长就有影响[3，18]。鉴于研究区为农作物种植耕地，农作物生长对盐渍化极为敏感，所以采用表 1分类方案进行盐渍化等级的划分，将土壤盐渍化分为5个等级。

表 1 盐渍化分级

2.2 光谱指数与特征波段选择

2.2.1 土壤光谱指数选择

土壤含盐量的变化影响着光谱的反射特征，总结这些特征可以反映土壤盐碱程度，因而可以通过选用前人建立的表征土壤盐渍化的光谱指数来反映盐渍化的发育情况。选取较为常用的6个直接和间接表征土壤盐分的光谱指数对顾屯耕地的盐分信息进行定量计算，将不同光谱指数的计算结果进行对比，利用效果最优的光谱指数反演研究区盐分含量，光谱指数的具体内容如表 2所示。

表 2 光谱指数

2.2.2 特征波段筛选

遥感影像不同波段之间反映的信息并不完全独立，不同波段传递的信息存在一定的相关性，相关性较大容易造成模型失真和不稳定[10]。同时不同波段对同一地物特征的描述程度(即波段的信息量)也是不同的，因此通过选取波段之间相关性小、单个波段信息含量大且对盐分含量较为敏感的波段组合理论上能够更好地反演盐分信息。首先提取高分一号影像中对应采样点的4个波段的波段值，与实验得到的盐分数据分别做相关性(R)分析[式(1)]。再对4个波段两两之间进行相关性分析，同时计算各个波段的方差(D)[式(2)]，数据方差大小反映了数据的离散程度，也即信息量的大小。根据计算结果筛选信息量较大、相互之间相关性较低且与盐分信息高相关的波段及波段组合即作为盐分反演的建模因子。

(1)

(2)

2.3 反演模型及判据

2.3.1 支持向量机

SVM是根据统计学习理论，以结构风险最小化原则为理论基础从线性可分扩展到线性不可分的一种机器学习方法，通过核函数来实现非线性映射来处理非线性问题，对于小样本情况下能取得较好的效果，能够较好地解决非线性回归问题。

设定研究数据集Q={(xi，yi)，i=1，2，…，n}，其中，xi为输入的波段变量组合数据，yi为训练样本的实测盐分数据，以此建立的SVM模型为

(3)

式(3)中：K(xi，x)为内积核函数。 SVM共有5种核函数类型，选择高斯径向基核函数

(4)

式(4)中：P为核函数的核宽度。

SVM学习的准确性取决于参数的选择，主要参数为惩罚参数C、不敏感带损失函数g和核函数的核宽度P[23]，其中g=1/(2P2)。参数C与g分别控制模型的泛化性能与回归误差，为确保模型精度，基于训练样本集交叉验证的网格搜索法进行参数寻优，由最小均方差的原则确定参数C与g。

2.3.2 判据评定

利用14个验证样本对构建的反演模型进行精度验证。分别计算各模型反演值与验证样本的均方根误差(root mean square error，RMSE)[式(5)]、决定系数R2与相对分析误差(relative prediction deviation，RPD)[式(6)]，当R2越接近于1、RMSE越小、RPD相对较大时，表明模型的预测效果较好。

(5)

(6)

3 结果与分析

3.1 基于盐分指数的盐渍化反演模型

利用高分一号影像在ENVI 5.3软件bandmath功能中将所选的6个光谱指数进行波段计算并进行归一化处理，在ARCGIS 10.6软件中提取计算结果值，与筛选的建模样本定量拟合(图 2)，可以看出6种光谱指数与实际的盐分含量都有一定的相关性，其中S3与盐分含量正相关[图 2(e)]，R2达到0.716，拟合度最好；S1与S2与盐分含量负相关，R2分别为0.680、0.632[图 2(a)、图 2(c)]； DVI的拟合度为0.507[图 2(f)]；NDVI与SI2与盐分含量的相关性不高[图 2(d)、图 2(b)]，R2只有0.156和0.061，无法用来拟合盐分含量。基于定量分析结果，选择S1、S2、S33个拟合度较高的指数建立光谱指数反演模型，模型具体如表 3所示。利用3个反演模型对上一步光谱指数计算结果进行计算，提取反演盐分值与验证样本进行反演精度评估[式(5)、式(6)]，结果如表 3所示。S1与S2模型的验证集R2分别为0.587与0.456，反演效果一般，难以用来进行研究区的盐渍化分布评价。而S3模型R2达到0.742，同时RMSE较小、RPD最高，表明S3具有较好的反演效果，S3模型能在一定程度上反映研究区的盐渍化特征。

图 2 指数与盐分定量关系Fig.2 Relationship of index and salt quantity

图 3 SVM模型评价精度Fig.3 Accuracy of SVM model evaluation

3.2 基于特征波段融合的盐渍化反演模型

高分一号影像共有4个多光谱波段，对土壤采样点对应的4个波段数值及盐分含量进行两两相关性分析[式(1)]得到相关性矩阵，如表 4所示。波段1与波段2相关性最为显著，相关系数达到0.690，波段1与波段4、波段3与波段4、波段1与波段3相关性较低，表明两两之间具有一定的独立性，传达的信息冗余程度较低。其中波段1与波段4相关性最低，且负相关，为-0.200。

表 3 反演模型及精度参数

表 4 相关性矩阵

表 5 信息量统计表

表 6 最优参数

表 7 SVM反演模型

按式(2)对各波段进行信息量统计，结果如表 5所示。波段1与波段4的方差最大，波段2方差最小，说明波段1与波段4所包含的信息量相对而言更大。从表 4可知，单波段与盐分含量的相关性并不特别显著，最大仅为0.413。综合因子需包含足够信息量且相关性较小的原则，选取两两之间相对独立的波段1与波段3、波段1与波段4及波段3与波段4作为评价顾屯流域耕地盐渍化的3个自变量组合。将3个模型变量归一化后结合筛选的42个训练样本在MATLAB平台利用交叉验证的网格搜索法分别进行参数寻优，寻优结果(表 6)即代表在相应参数下各变量组合在SVM模型中能够最大程度地描述出数据特征，建立最优的非线性回归模型。基于所得最优参数分别构建3个波段组合的反演模型[式(3)、式(4)]，对研究区进行盐分空间分布的反演，提取反演结果进行精度验证[式(5)、式(6)]，结果如表 7与图 3所示。P1P4、P3P4组合训练集R2分别为0.741、0.725，而预测R2仅0.679、0.672，表明对于耕地盐渍化反演效果略差。P1P3组合建模集R2高达0.901，预测R2也达到0.8674，RMSE为0.1052，RPD为2.931，各项指标表明该模型盐分反演的效果最好。

图 4 盐渍化分布Fig.4 Distribution of salinization

3.3 两种模型反演结果对比分析

从评价精度上，由盐分指数S3构成的反演模型，R2达到0.742；由P1P3组合构建的SVM模型在所有反演模型中精度最高，R2为0.867，相较S3模型具有更优的预测能力，能够高精度反演耕地的盐分分布情况。这说明波段1与波段3是最优的波段组合，同时SVM能够很好地挖掘数据特征进行分析建模。利用S3反演模型(表 3)和P1P3支持向量机反演模型(表 7)分别对顾屯流域耕地区盐分进行空间分布反演，并按照表 1分级规则对盐分进行分级，获得耕地盐分空间分布情况，如图 4所示。

表 8 反演结果统计表

对两个模型评价结果进行统计，结果如表 8所示，从S3模型的评价结果看，顾屯沟内约70%耕地都受到一定程度的盐渍化影响，其中轻度盐渍化占比最大，占了顾屯主沟56%，中度盐渍化土壤为17.39%，重度盐渍土仅有0.03%。SVM模型得出的结果显示，轻度盐渍化占43.52%，中度盐渍土占29.82%，重度盐渍土占0.04%。对比可以发现，两个模型反演所得正常土壤与重度盐渍化土壤的总量相差不大，总量上的主要差别在于SVM模型得出中度盐渍土要多大概10%，相应的轻度盐渍化相对要少10%，即SVM模型表明研究区的盐渍化趋势更强，从图 4也可以看出SVM模型中度盐渍土占比较多。两种模型的反演结果(图 4)均显示，在沟内上游区域与支沟内盐分含量明显高于下游，这可能与沟内高程与水文条件等因素有关。盐分空间分布结果表明研究区经工程后的耕地在一定程度上受到盐渍化的影响，盐渍化会对耕地农作物的生长与工程后土地的利用效益产生负作用。

4 结论

土壤盐渍化问题影响耕地的种植和治沟造地工程的效益。通过高分一号影像和室内实验数据对顾屯流域耕地的盐渍化反演，初步获得以下结论。

(1)选取的6个表征盐分光谱指数与实测盐分信息具有一定的相关性，其中S1、S2、S33个指数与实际盐分成分有较好的相关性。通过定量反演，S3模型也表现出了较好的反演精度，常用的盐分指数对于研究区盐渍化评估有一定效果。

(2)通过信息量、相关性等特征筛选出高分一号影像P1P3、P1P4与P3P43组特征波段组合，发现P1P3可作为盐渍化反演最优波段组合。

(3)利用SVM的机器学习方法对耕地盐分含量进行反演，模型精度达 0.867。该模型可为黄土耕地盐渍化的监测提供指导。

(4)顾屯流域70%的耕地正受到盐渍化影响，重度盐渍土极少，中度与轻度盐渍土成分较高，分别达到43.52%与29.82%。这说明研究区沟内耕地的盐渍化发育态势较为显著，盐渍化对于耕地农作物的种植有害。当地的劳动者可按照盐渍化的分布情况合理的选择农作物的种植区域与种类，对于盐渍化程度相对较高的区域应该及时治理，同时对于轻度盐渍化区域也应及时抑制土壤盐分含量的增加。

需要指出的是虽然对耕地盐分进行了较高精度的反演，但也存在一些缺陷：高分一号影像光谱分辨率有限，仅有的4波段对于更准确地描述盐分信息来说较为困难；仅研究了顾屯耕地的盐渍化空间分布情况，对于盐渍化的形成机制及工程后效应没有做深入的讨论，要解决治沟造地工程带来耕地盐渍灾害问题，还一定要弄清盐渍化的形成与工程活动、地形地貌、水文的关系，这是今后的研究方向。