航空高光谱反演耕地土壤重金属分布特征——以苏北灌河地区为例

徐明钻,梁 森,石剑龙,季 岩,黄 岩,梁胜跃,严维兵

(1.江苏省地质勘查技术院,江苏 南京 210049;2.江苏省地质矿产勘查局航空对地探测技术研究中心,江苏 南京 210049)

21世纪以来,随着中国人口的快速增长,各类产业蓬勃发展,对工业产品的不断需求及日益增大的能源消耗,导致大量重金属和有机污染物不断向各种环境介质排放并积累,排污能力已超过土壤、水体和大气系统的自净能力[1],因矿山选冶、化工排污、农药化肥不合理使用等造成的耕地土壤重金属污染问题尤其突出[2-3]。由于重金属具有持久性和不可降解性,极易对生态环境造成污染,严重影响动植物生长发育和农产品质量安全,甚至危害人类健康[4-5]。目前,耕地土壤重金属污染问题已经成为社会关注的热点和焦点,受到各级政府的高度重视,因此,迫切需要开展耕地土壤重金属污染调查与修复治理工作。

航空高光谱遥感技术广泛应用在监测土壤环境、评价土地质量、编制土壤类型图以及定性半定量预测土壤成分等研究中,利用高光谱遥感技术获得的土壤表层信息及其空间分布信息,可达到探测土壤成分、预测土壤理化性质的目的[6]。目前,国内外众多学者进行了反射光谱反演土壤重金属元素含量研究,例如KOOISTRA L等[7]利用河滩土壤可见光-近红外光谱和土壤重金属Zn、Cd含量建立PLSR反演模型,定量分析了河滩土壤成分和重金属含量;WANG J J等[8]利用室内测定光谱对宜兴农田土壤重金属Pb、Zn、Cu建模反演,获得了显著效果。尽管高光谱反演土壤重金属含量方面已取得了一定成果,但这些成果大多数是基于室内光谱或小区域范围进行的研究工作。本文探索机载航空高光谱在大比例尺、大面积尺度下快速反演耕地土壤重金属含量和空间分布的技术方法及工作流程。以苏北灌河地区为例,利用机载航空高光谱遥感方法快速完成了研究区调查工作,通过土壤重金属在可见光-近红外波段光谱的吸收特征,分析土壤光谱特征与污染元素的定量关系,建立土壤污染元素含量估算模型,进一步评价土壤污染程度,并与传统地球化学方法进行对比分析,验证评价结果的可靠性,为该区农业生产及生态环境保护提供参考。

1 研究区概况

研究区位于苏北灌南县与响水县交界处,濒临黄海,新沂河、灌河由西南向东北流经该区入海,面积为582 km2。根据地层成因类型、岩石组合及其地貌分布等特征划分, 该区为苏北冲海积平原的组成部分,大部分属于沿海平原区,靠近响水县城部分属于徐连丘岗平原区,无基岩出露,松散层厚度近千米。区内土壤普遍含较高的可溶性盐分,但受成陆先后、海水侵袭、黄河决口等影响,土壤发生了复杂演变,同时经过多年治水、垦殖、耕作、改良等人工改造,区内土壤大致可以分为3个类型:① 外滩砂土、两合土区;② 里滩砂性盐土、壤性盐土区;③ 圩里油泥土区。灌河沿东部浅部土体以砂性土为主,透水性较强,抗污染能力弱-中等;其余地段浅部土体黏性土为主,且厚度较大,透水性弱-微,抗污染能力强。

灌河东岸有陈家港化工园,西岸有连云港堆沟化工园区、临港产业区,分布数十家高污染化工企业,周围大气、土壤及灌河均受不同程度污染,农业面源污染也较严重。田爱军等[9]研究表明,灌河流域周边农田化肥平均使用量达52 500 kg/km2,高于全国农田化肥平均使用量(21 855 kg/km2),加大了区内土壤污染负荷,同时该区耕地主要就近采用地表水进行灌溉,加剧了耕地土壤污染程度。研究区西南部、灌河以北耕地土壤重金属元素含量明显高于其他地区[10],因此,本文将研究区划分为一般区和重点区,具体位置见图1。

图1 研究区卫星影像图

2 研究方法

2.1 高光谱成像仪及数据获取

本次工作采用自主集成的机载航空高光谱遥感测量系统搭载在运-12固定翼飞机进行数据采集。该系统主要由中国地质调查局南京地质调查中心研发的国产成像光谱仪HSM400A和HSM1000A、三轴稳定平台GSM4000、惯性导航系统KVH1750、数据收录控制系统和导航控制系统。可获取亚米级空间分辨率影像、波长范围400～1 000 nm和1 000～2 500 nm高光谱数据,总波段数543个。HSM400A和HSM1000A主要技术指标如表1所示。

表1 国产成像光谱仪参数

2017年11—12月,在庄稼收割后翻土期开展了野外数据采集工作,完成了研究区582 km2的野外航空高光谱数据信息采集,获取的HSM400A(400～1 000 nm)光谱分辨率为2.42 nm,空间分辨率为1.0 m;HSM1 000A(1 000～2 500 nm)光谱分辨率为6.98 nm,空间分辨率为1.75 m。地面地物光谱采样161条,地面土壤样及分析测试128个。

2.2 数据预处理

机载航空高光谱影像具有较高的光谱分辨率和空间分辨率,机载平台下的数据采集面临复杂的外界环境变化,光谱存在条带噪声、几何畸变、信息冗余等问题,针对机载航空高光谱数据的预处理尤为重要。飞行前,利用单色仪和积分球对成像光谱仪进行实验室光谱标定和辐射标定,建立辐射定标模型;获得辐亮度数据后,使用基于辐射传输理论的MODTRAN模型进行像素级大气校正,消除水汽、二氧化碳、臭氧等分子和气溶胶散射效应[11-12],矫正目标像元和邻近像元的交叉辐射效应,获取地表反射率数据;利用HyperPic 软件进行几何校正计算。由于HSM400A和HSM1 000A数据的空间分辨率不同,需要对HSM400A数据进行空间分辨率重采样至1.75 m,将两种数据进行像元几何校正对齐,再对HSM400A和HSM1 000A光谱进行波段匹配叠加,生成波谱范围为400～2 500 nm的高光谱数据。获取的影像光谱反射率曲线见图2。

图2 影像光谱反射率曲线

2.3 提取特征波段算法

与多光谱相比,高光谱数据波段较多,可以提供更丰富的细节信息,但波段之间具有较高的冗余性,冗余的波段参与模型计算时会影响模型的稳定性和反演精度,因此,需要筛选出与重金属含量变量相关性较高的特征波段。本文采用竞争性自适应加权采样方法(Competitive Adaptive Reweighted Sampling,CARS)筛选土壤光谱波长变量,选择与土壤实验室元素分析结果相关的波长变量[13-14],提取的最大因子数由蒙特卡罗(MC)法交叉验证决定[15],MC采样次数为400次,选择的波长变量子集采用10折偏最小二乘法交叉验证建模,根据最小均方根误差值确定最优的波长变量子集。

2.4 建模分析与精度评价

利用选择的特征波段与重金属元素化验分析数据以BP神经网络(Back Propagation Neural Networks,BPNN)回归模型进行建模反演,利用交叉验证方法验证模型的预测能力。由于BP神经网络模型不需要具体的解释模型和理论推导,不需要具体的计算公式、繁杂的参数选取和太多的专家经验,同时具有很强的非线性处理能力和自适应特点,可以满足建立光谱特征与污染元素之间的对应关系。BP神经网络模型是一种多层前向型反馈隐含层与输出层构成,将偏差通过隐含层向输入层逐层逆向传递,以此逐级修正偏差[16-17]。本文通过特征光谱和元素浓度建立相关性模型,模型采用6层线性隐含层,其神经元分别为20、16、10、6、4、2,该项参数需根据反演结果进行调整,优化函数采用随机梯度下降法,通过设定迭代次数或目标精度终止循环得到最终模型[18-20],模型结构见图3。

图3 BP神经网络(BPNN)结构示意图

研究区反演模型共选择128个样本,选择方式为采样点去除异常点后,根据浓度排序将样本以2∶1划分为训练集和验证集,进行建模和验证。本文采用决定系数R2、均方根误差RMSE、平均绝对误差MAE等3个参数来评价所建模型的优劣。R2值越接近1,说明建模精度越高,RMSE值和MAE值越小,说明预测能力越强[21]。

3 结果与分析

3.1 特征波段提取

特征波段筛选结果见表2。首先,将土壤原始光谱进行降维处理并剔除掉噪声较大的波段,并采用Savitzky-Golay滤波进行平滑去噪处理,相比于其他滤波器,该滤波器最大的特点在于滤除噪声的同时可以确保光谱数据的形状及宽度不变[22]。然后,采用CARS选择法对土壤光谱波长变量进行筛选(图4),选出均方根误差(RMSE)最小子集为最优变量子集。设定蒙特卡罗采样为400次,由交叉验证得到的RMSE趋势变化可知,随着采样次数的增加,RMSE值由大到小,再逐渐增大。当Cd、Zn、Cu、As采样次数循环到246次、216次、251次、300次时,RMSE值最小,所选的波长子集最优。由表2可知,各重金属元素最优的波长子集包含的波段数为12～15个,波段主要集中在400～900 nm和1 300～2 300 nm,涵盖了土壤重金属特征吸收区,与袁自然等[23-25]研究的特征波段范围较相近。

表2 特征波段筛选结果

图4 自适应加权采样方法(CARS)选择变量-均方根误差变化图

3.2 模型训练与精度评估

耕地土壤重金属含量建模精度结果如表3所示。重金属As、Cu、Zn、Cd训练集RC2为0.915、0.893、0.874、0.911,均达到了0.8以上;验证集RP2分别为0.625、0.637、0.603、0.610,RMSE和MAE均在较小范围内,预测结果偏差较小。验证集RP2为0.603～0.637,低于袁自然等[23-25]获得的土壤重金属含量预测的相关系数,主要由于本文研究采用的高光谱数据是基于机载平台获取,其复杂的外界环境变化导致光谱质量有所损失。

表3 耕地土壤重金属含量建模精度结果

结合BPNN模型反演的预测值与实测值散点图(图5)可知,训练集和验证集的样本实测值和估算值均在1∶1线附近,特征光谱与土壤重金属含量之间具有较高的相关性,且4种重金属元素模型均具有相对较高的预测R2和相对较低的RMSE和MAE,说明模型具有较高的稳定性和预测精度。

图5 BP神经网络(BPNN)模型反演预测值与实测值散点图

3.3 与地球化学调查成果对比分析

参照相关技术规范[26-28],土壤污染元素浓度高光谱反演估算结果可划分为背景区、低风险区、高风险区,从宏观趋势上评价土壤污染状况。参照土壤地球化学背景值计算方法,根据高光谱反演污染元素含量值,利用多次迭代手段确定高光谱反演结果背景值。高风险值为规范中农用地土壤污染风险筛选值,低于背景值的为背景区,大于背景值小于高风险值的为低风险区,大于高风险值的为高风险区。反演结果划分区域限值如表4所示。

表4 反演结果划分区域限值

根据表4对研究区重金属污染光谱反演数据进行分级划分,结果显示研究区主要污染元素为As,且集中在重点区;Zn分布主要为背景值区和低风险区;Cd、Cu存在局部高风险区,由于受地物或植被影响,Cd、Cu个别高风险点未反映。针对重点区As高光谱反演空间分布图(图6)做进一步评价,并与1∶1万地球化学调查As浓度空间分布图(图7)做对比分析,从土壤重金属As分布趋势上验证光谱反演结果的可靠性。

图6 As高光谱反演空间分布图

图7 1∶1万地球化学调查As浓度空间分布图

由图7可知,重点区绝大部分地段As含量为10～15 mg/kg,>26 mg/kg的高值区主要呈NE向-SW向条带状分布在研究区中部和南部。结合高光谱反演结果(图6),研究区中部和南部NE向-SW向A区、B区和C区反演结果显示为高风险区,其范围、位置和走向与图7中As实测异常相对吻合,另外实测主异常带西南侧还存在NW向As 弱异常带,而相应的高光谱反演解译结果也显示出对应的NW向As弱异常带。高光谱反演结果在趋势上与传统地球化学调查结果吻合度较高,验证了高光谱反演结果的可靠性。

4 结论

(1)航空高光谱遥感技术在耕地土壤重金属污染快速调查工作中,能较好地评价耕地土壤重金属污染程度,建立的航空高光谱遥感土壤和水体污染快速调查技术方法和工作流程具有一定可行性,为探索大尺度应用航空高光谱遥感方法调查耕地土壤重金属污染状况奠定了基础。

(2)基于机载航空高光谱数据在获取时,光谱质量有所损失,模型预测的决定系数虽然略低于实验室数据,但也显示了预测值和实测值之间存在较为显著的相关性,说明建立的模型具有相对较好的回归性能。

(3)相对于传统地球化学方法调查耕地土壤获得以点带面的数据成果,高光谱反演解译结果显示的耕地土壤As高风险区异常更全面精细,弥补了传统地球化学方法离散采样点间的空白区。航空高光谱在大面积的耕地土壤重金属调查工作中,可以在一定程度上实现宏观、实时、快速、有效的调查。