香蒲砷含量高光谱估算模型研究

余哲修 张 超 任 伟 张 一 牛晓花 王 妍 刘云根

( 1. 北京林业大学省部共建森林培育与保护教育部重点实验室，北京 100083；2. 西南林业大学林学院，云南 昆明 650233；

3. 西南林业大学生态与水土保持学院，云南 昆明 650233)

20 世纪90年代以来，我国湿地重金属污染问题加剧，不仅对湿地生态系统造成了严重破坏，对周边植被生长亦产生深远影响[1]。利用高光谱技术定量反演分析湿地植物的理化参数时序变化特征，能有效避免传统湿地植被调查中成本高、耗时费力等弊端，为今后精确估算、定量评价、有效保护及合理利用湿地资源提供了技术手段和数据支撑[2]。

目前，应用高光谱技术探究农作物、土壤和非湿地环境植物在重金属污染胁迫下的光谱特征变化已成国内外为研究的热点和前沿[3-4]，而在湿地方面的研究报道较少。光谱分析技术能够准确描述和定量表达重金属污染胁迫下湿地植物的生理生化特征变化，可发挥高光谱遥感技术的高效性、经济性和准确性等优势[5]。陈思宁等[6]通过锌添加土培白菜（Brassica campestris）盆栽实验研究了锌胁迫的白菜叶片光谱特征，所提取的3 种特征光谱因子能够评价白菜锌污染胁迫的程度。Slonecker 等[7]在不同砷浓度梯度下栽培蜘蛛人凤尾蕨（Pteris cretica mayii）、井边栏草（Pteris multifida）和波士顿蕨（Nephrolepis exaltata），基于其反射率一阶导数构建的砷浓度估测模型较为合理准确。刘美玲等[8]利用ASD 光谱仪采集了不同铅污染水平的水稻（Oryza sativa）反射光谱，通过小波变换、分形分析和模糊数学的结合实现了利用反射光谱估算水稻重金属污染胁迫状况。Li 等[9]基于Hyperion 影像和地面实测光谱数据，结合建立的高光谱变量，能够较为准确地监测研究区中砷胁迫下的水稻冠层叶绿素含量。屈永华等[10]利用HJ-1A HSI 影像研究了江西省德兴铜矿分布区周边植被对铜胁迫的光谱响应，其研究结果可为矿区铜污染的预防、治理和监测提供技术支持。Gao 等[11]通过高光谱数据和支持向量机回归对太湖湿地土壤的锌、砷和铜含量进行估测，该方法对均匀湿地土壤中重金属含量估算效果较好。

湿地作为湖泊与陆地之间重要的过渡带，通过湿地植物重金属含量的研究，可以为湿地生态系统污染治理和修复提供思路，定量拟合植被反射光谱与主要生理生化指标之间的关联关系，可达到宏观监测湿地植被生长状态的目标。众多学者研究的植物多数为易受重金属影响其生长的植物，香蒲（Typha orientalis）作为一种典型湿地植物，并且对重金属具有富集能力，在修复生态系统方面受到重视[12-13]。针对其在重金属胁迫环境下反射光谱特征变化及响应的研究较少，故本研究基于不同砷含量胁迫条件下香蒲的反射光谱特征，筛选出与砷相关性高的光谱变量，对香蒲砷含量进行高光谱估算建模，旨在为今后利用高光谱遥感技术和地面高光谱测定协同精确反演湿地植物的重金属含量提供参考。

1 材料与方法

1.1 实验设计

2017年4 月在昆明郊区苗圃采购的当年生香蒲幼苗为实验材料，根据《土壤环境质量 农用地土壤污染风险管控标准（试行）》（GB 15618—2018）[14]、云南土壤砷含量背景值和2015年10月在阳宗海南岸湖滨湿地采样的结果[15]，实验采取盆栽土培的方式，设置9个砷胁迫梯度，土壤全部采用当地的红壤，粉碎块状土块，混合均匀，熟化2 周。

实验于2017年5月进行，砷浓度梯度分以添加Na2HAsO4·7H2O 固体的形式设置。每个梯度设6个重复。黑色培养桶（高32 cm，内径44 cm，外径48 cm，底径32 cm）装入10 kg 处理后的土壤，土壤厚度约15 cm，加入10 L 水，积水深度10 ～ 15 cm，模拟野外湿地底泥积水环境。搅拌均匀后对底泥进行污染，期间间隔一定时间对其重复搅拌，保证底泥能够充分污染。每个重复6 株香蒲，另设实验对照1个，共10个梯度，具体设置如表1 所示，其中梯度编号规则为：胁迫类型-植物类型-胁迫浓度。实验采用统一的水肥管理，保证幼苗生长健壮、无病虫害、长势基本一致。将混匀后的土壤装入培养桶中，在各培养桶外壁分别标识表1 所示的砷浓度梯度编号。光谱采集后将香蒲取样送至实验室，使用原子荧光法测定香蒲叶片的砷（As）含量。

表1 不同砷浓度梯度下样本设置Table 1 Number and concentration of samples in different arsenic gradients

在室外采集光谱的同时，使用CCM-200 Plus叶绿素仪同步测定香蒲叶片叶绿素含量指数（CCI），该仪器通过测量叶片在653 nm 和931 nm处叶绿素不同吸收率确定CCI 值。室内测定砷含量的同时，测定了叶片鲜质量、含水率和磷含量，鲜质量通过天平测定，含水率根据鲜质量和烘干后得到的干质量计算得到，磷含量测定采用钼锑钪比色法。

1.2 光谱采集与处理

利用ASD Field Spec 3 野外光谱仪采集光谱数据（波长为350～2 500 nm）。测量时在晴朗、无风、少云条件下进行，测量时间为10：00—15：00。光谱采样时，将香蒲叶片展平，仪器探头垂直于叶片，每个样本测量30 次光谱，视实际环境随时对仪器定标。使用ViewSpec Pro 软件对原始光谱数据进行处理。受光照强度、风、空气湿度等影响，需要对原始光谱数据进行筛选，剔除异常数据，对各浓度梯度的光谱做平均化处理和平滑处理，使光谱曲线更准确地表现香蒲在各砷浓度梯度下的光谱特征。在整段光谱曲线中，剔除被噪声影响的1 350～1 400、1 750～1950 nm 两部分波段。

1.3 光谱变化分析

1.3.1 光谱导数

对光谱曲线作导数变换，导数变换不仅能凸显不同光谱曲线之间的差异，亦能消除部分大气、土壤等影响，更准确反映植被光谱特征，更好地提取光谱特征波段[16-18]。根据光谱反射率获得一阶导数和二阶导数光谱，公式如下：

式中： FDRλi为波段i和i+1 波长值间的一阶微分反射比； SDRλi为波段i和i+2 波长值间的二阶微分反射比；Rλi、Rλi+1、Rλi+2分别为波段i、i+1、i+2处的光谱反射率； Δλ为波段i至i+1 的波长值。

1.3.2 “三边”参数

三边参数相比原始光谱可以更准确反映植被的光谱特性[19-21]，从原始光谱和光谱一阶微分中提取20个“三边参数”进行分析，如表2 所示。

表2 “三边”参数变量Table 2 Variables of 'three edges' parameters

1.3.3 窄波段植被指数

植被指数通过多个波段进行代数运算，在提取植被信息的灵敏度方面相比单个波段较高[22]，常用的植被指数主要有差值植被指数（DVI）、归一化植被指数（NDVI）、比值植被指数 （RVI），在原始计算方式中参与计算的波段为红光波段和近红外波段，窄波段植被指数突破了原始计算方式波段的限制，通过全波段参与运算，使得植被指数的组合更多，同时能从中提取与植被理化参数更密切相关的波段组合。反射率经过一阶导数和二阶导数变换后，能够突出更多的植被光谱信息，因此将导数变换后的反射率值也参与运算，构建的窄波段植被指数见表3。和分别为第i和第j波段的反射率一阶微分值；和分别为第i和第j波段的反射率二阶微分值。

表3 窄波段植被指数变量Table 3 Variables of narrow band vegetation index

1.4 建模方法

1.4.1 光谱特征变量筛选

分别分析反射率、反射率一阶微分、“三边”参数和反射率二阶微分与砷的相关性，优选出与砷相关性高的变量用于估算模型建立。在60个香蒲样本中，随机选择45个用于建模，其余15个作为验证数据。

1.4.2 模型选择

单变量回归分析将显著性最高的光谱特征变量作为自变量用来反演理化参数，常用的单变量拟合模型如表4 所示。

多元回归分析将植被理化参数作为因变量，与理化参数相关性较高的波段、波段组合等其他光谱特征变量作为自变量来反演植被理化参数。以理化参数为因变量，光谱特征变量等作为自变量，构建多元逐步回归模型，如式（3）所示。在建立模型时对自变量逐个筛选，最终保留与理化参数相关性较高的自变量。

式中：y为估算的理化参数；an为第n个变量的系数；xn为第n个变量；ε为常数。

表4 单变量拟合模型Table 4 Single variable fitting model

1.4.3 模型检验

通过决定系数（R2）和均方根误差（RMSE）评价模型精度。模型决定系数是衡量建模效果优劣的常用指标，是衡量模型所推算的估算值与真实值间的差异化程度，R2在0 ～ 1 之间，值越大模型预测精度越高。均方根误差可以反映估算数值偏离实测值的程度，RMSE 值越小模型的精度越高。2 种评价指标计算方法如下：

式中：yi为第i个样本的实测值；为第i个样本的估算值；y¯为所有样本实测值的平均值。

2 结果与分析

2.1 光谱特征分析

2.1.1 原始光谱分析

不同砷浓度下香蒲的反射光谱曲线见图1。各浓度梯度下香蒲的反射光谱特征均保持了健康绿色植被反射光谱特征，从10 条光谱曲线中分段进行分析。由图1 可知，在可见光波段，A-T-2 在550 nm 的反射峰附近反射率最高，A-T-8 最低，在670 nm 的吸收谷附近，各浓度梯度下香蒲的反射率几乎相同；在近红外波段具有强烈反射特征，各梯度的香蒲在近红外波段有明显区别，A-T-8 的反射率最低，A-T-5 的反射率最高；在1 400～1 750 nm 的短波红外范围内，A-T-8 的反射率最低，1 950～2 450 nm 的短波红外范围内，AT-8 和A-T-9 的反射率较低，说明高浓度砷对香蒲的水分吸收产生一定影响。

图1 不同砷浓度胁迫下香蒲反射光谱特征Fig. 1 Spectral characteristics of T. orientalis in different concentrations of arsenic stress

2.1.2 一阶微分光谱分析

植被光谱特征主要体现在可见光-近红外波段，因此选取350～1 350 nm 对香蒲光谱反射率一阶微分值进行分析。由图2 可知，经过一阶微分变换后，各浓度梯度下香蒲的光谱反射特征的区别在红光和近红外波段体现出来。红光波段的峰值反映了原始光谱红光吸收谷变化最快的波长位置，各浓度梯度下香蒲的吸收谷有明显区别，AT-2 的一阶微分值最大，A-T-8 最小，表明A-T-2 叶绿素对红光的吸收最快，A-T-8 最慢，说明砷对香蒲叶绿素吸收红光有一定影响；由于植物叶片叶肉海绵组织的特性，在近红外波段使入射光形成多重反射，而950 nm 和1 150 nm 附近有水、氧气的窄吸收带[23]，因此形成了光谱曲线形态波状起伏的特点[24]，从图中可以看出950 nm 和1 150 nm 各曲线有明显区别，说明砷对香蒲的叶片组织产生了一定影响，使得香蒲对水分和氧气的吸收有明显区别。

图2 不同砷浓度胁迫下香蒲反射光谱一阶微分特征Fig. 2 First order derivative spectral of T. orientalis in different concentrations of arsenic stress

“三边”参数可以定量的对光谱分析，相比原始光谱更准确反映植被光谱特性。各梯度香蒲反射光谱的“三边”参数见表5。红边是反映植被生长状况的指示因子，植被长势较好时红边向近红外波段移动，植物受到胁迫或衰老时红边向蓝光方向移动[25]，由表5 可知，各梯度香蒲的蓝边、黄边和红边位置差异不明显，说明香蒲在砷胁迫条件下仍然能表现出很强的耐性；从绿峰和红谷上看，波长位置和反射率差异不明显，表明砷对香蒲的叶绿素没有影响，但在绿光和红光波段均表现出典型的反射和吸收特征；从面积变量和比值变量上，同一面积变量之间差异不显著，而比值变量略有差异，说明通过比值变量能够看出各浓度梯度砷胁迫下的香蒲反射光谱间的差异。

2.1.3 二阶微分光谱分析

从一阶微分曲线上看，各梯度香蒲的光谱特征差异较小，二阶微分相比一阶微分可以更有效地减弱大气、背景噪声等因素，从而分辨重叠光谱，二阶微分曲线的凹凸形态可以表现出原始光谱的反射和吸收特征，凹形代表波段在吸收区域，凸形代表波段在反射区域[26-28]。以350～1 350 nm 为主要分析波段，由图3 可知，经过二阶微分变换后原始光谱的反射和吸收特征相比一阶微分更好的表现出来，可见光波段中510、550、690 nm 和730 nm 为主要的特征波段，在近红外波段中900～950 nm 和1 100～1 150 nm 有明显区别，更好的体现了原始光谱近红外波段对水和氧气窄吸收带的光谱特征。

表5 “三边”参数统计Table 5 The statistics of 'three edges' parameter

图3 不同砷浓度胁迫下香蒲反射光谱二阶微分特征Fig. 3 Second order derivative spectral of T. orientalis in different concentrations of arsenic stress

2.1.4 理化参数统计分析

通过光谱一阶导数、二阶导数变换和“三边”参数分析了不同浓度砷胁迫下的香蒲的光谱特征，但光谱特征差别较小，说明从反射光谱上表现出了香蒲对砷的富集能力，结合其他的生理生化参数有助于分析其光谱响应。理化参数均值和标准差统计如表6 所示。

由表6 可知，CCI、磷含量和含水率的均值较为稳定，鲜质量的均值和标准差波动较大，说明香蒲叶片的磷含量和鲜质量对不同的砷含量响应有所不同，而理化参数之间互相的作用对其光谱响应也有影响，因此需要进一步分析理化参数间的相关性，如表7 所示。由表7 可知，含水率与磷和鲜质量达到极显著相关（P＜0.01），香蒲作为一种挺水植物，水环境对其生长有显著影响，磷则影响其生长与品质，对其他营养元素的吸收有一定影响，叶绿素与植物光合作用相关，在理化参数的综合作用下影响其鲜质量。砷与其他4个参数的相关性较弱，均未达到显著相关（P＜0.05），表明砷对香蒲的叶绿素、磷、鲜质量、含水率影响较小，进一步说明了香蒲对砷有较强的富集能力。

表6 理化参数均值和标准差统计Table 6 The mean and standard deviation of physicochemical parameters

表7 理化参数相关性分析Table 7 The correlation analysis of physicochemical parameters

2.2 砷含量估算模型构建与评价

2.2.1 相关性分析

分别将原始光谱反射率、一阶微分、二阶微分和“三边”参数与砷含量进行相关性分析。由图4 可知，原始光谱反射率所有波段与砷的相关性均为负相关，1 475个波段与砷的极显著相关（P＜0.01），在1 115 nm 相关性最高，相关系数为-0.443。原始光谱经过导数变换后，特征波段更加明显，一阶微分中290个波段达到极显著正相关（P＜0.01），947 nm 相关性最高，相关系数为0.435，425个波段达到极显著负相关（P＜0.01），713 nm 相关性最高，相关系数为-0.444。二阶微分中129个波段达到极显著正相关（P＜0.01），733 nm 相关性最高，相关系数为0.456，96个波段达到极显著负相关（P＜0.01），973 nm 相关性最高，相关系数为-0.494。

图4 反射率及反射率微分值与砷含量的相关性Fig. 4 The correlation between arsenic content and the value of reflectance or reflectance derivative

由图5 可知，在20个“三边”参数中有10个参数与砷含量显著相关（P＜0.05），其中3个参数达到极显著相关（P＜0.01），其中，SDr相关性最高，相关系数为-0.447。

图5 “三边”参数与砷含量的相关性Fig. 5 The correlation between arsenic content and'three edges' parameters

高光谱数据中的相邻波段之间具有高度的相关性，且相邻波段间提供的信息较为相似，造成了数据冗余，构建窄波段植被指数时，将光谱数据重采样为5 nm，再通过重采样后的光谱数据进行组合，一方面可减少信息重复，另一方面可减少运算量。重采样后共有373个波段，每种植被指数各有138 756个组合。

9个窄波段植被指数与砷含量的相关性见图6。从矩阵图形上看，基于导数光谱的相关系数矩阵离散程度较高，正负相关交替的范围较多，整体上看，红光波段、近红外波段和短波红外波段的组合相关性较高。从相关性上分析，由原始光谱任意两波段组合的NDVI 和RVI 相比原始光谱单波段与砷含量的相关性明显提高，说明窄波段植被指数组合能够明显提高香蒲对砷的响应。导数变换后组合的NDVI 和RVI 在相关性上又有明显提升，二阶导数相比一阶导数的组合相关性略有下降，但DVI、FDVI 和SDVI 最高相关系数绝对值均未超过0.5，说明差值指数的组合对砷的敏感性较低，FNDVI 和FRVI 最高相关系数均大于0.7，最佳组合为FRVI (1 065， 1 260)，相关系数为0.760，原始光谱经过导数变换后，能够消除部分背景噪声，因此FNDVI、FRVI、SNDVI 和SRVI 基于导数光谱构建的指数与砷含量表现出较好的相关性。

图6 窄波段植被指数与砷含量的相关系数矩阵Fig. 6 The correlation coefficient matrices between arsenic content and narrow band vegetation index

2.2.2 模型评价

将相关性最高的变量用于建立估算模型，估测模型评价指标见表8。基于FRVI (1 065， 1 260)建立的三次函数模型效果最佳，R2最高，RMSE最小，说明基于FRVI 的模型可用于砷含量的估算。基于SDVI (970， 925) 构建的估算模型效果较差，RMSE 最大，模型拟合精度最低。

对比分析普通导数变换与窄波段植被指数构建的估算模型可知，除DVI、FDVI 和SDVI 外，其他基于导数光谱变换结合窄波段植被指数的估测精度高于普通导数变换，有效地提升了香蒲砷含量的估测精度，R2最高可提升0.472。宽波段的RVI 和NDVI 由近红外和红光波段组合，近红外波段对植被的敏感性较高，而高光谱任意选取的两个波段的组合不受波段范围限制，获取的信息更加丰富，充分挖掘了植被的光谱信息，提高了光谱信息可利用程度，而导数变换有利于获取植被光谱细节信息，因此基于导数光谱构建的窄波段植被指数可以有效提高普通导数变换对香蒲砷含量的估测精度与稳定性。结合普通光谱变换与窄波段植被指数有效地挖掘了植被光谱内的可用信息，削弱了大部分背景噪声的影响，同时提高了光谱对砷含量的敏感性，可以为利用高光谱遥感反演植被理化参数提供参考。

表8 香蒲砷含量估测模型Table 8 Estimation models for arsenic content of T. orientalis

3 结论与讨论

本研究对不同砷浓度梯度环境下培养的香蒲反射光谱进行分析，并对其砷含量进行估算建模，结论如下：

1）光谱特征方面，不同砷浓度梯度生长的香蒲均表现出健康绿色植被的反射光谱特征，但第8 和第9 梯度的香蒲在1 950～2 450 nm 范围内在反射率低于其他梯度，说明高浓度砷对香蒲的水分吸收产生一定影响；通过一阶导数和二阶导数变换后，相比原始光谱更能凸显各梯度的香蒲反射光谱在近红外波段间的差异，“三边”参数的比值变量可以更好地定量分析光谱差异。

2）原始反射率与砷含量间全部为负相关，反射率一阶微分、“三边”参数和反射率二阶微分与砷含量间均表现为正负相关波动的变化，相比原始光谱可以更好地提取特征波段，提高反射率与砷含量的相关性，经过二阶导数变换后，973 nm的反射率二阶微分与砷含量的相关性最高，相关系数为-0.494；构建的9 种窄波段植被指数能有效提高反射率与砷含量的相关性，相关性最高的指数组合为FRVI (1 065， 1 260)，相关系数为0.760。

3）分别通过反射率、反射率一阶微分、“三边”参数、反射率二阶微分和植被指数建立砷含量的估算模型，综合比较得出，基于FRVI (1 065，1 260) 构建的回归模型估算效果最佳，决定系数为0.793，均方根误差107.457。

通过原始光谱曲线和变换后的光谱曲线特征能够较为准确地描述香蒲在不同砷浓度环境下的差异，在以后的研究中可以尝试利用室内光谱设备采集光谱数据，避免水汽、太阳光照变化的影响，并且环境时间可控，能获得更加准确的光谱数据。在未来的实验设计上，可以从香蒲多个生长期的理化参数和反射光谱，或通过野外实验，进一步分析不同生长期光谱和砷含量的相关机理。在光谱分析方面，仅通过原始光谱反射率、光谱导数、植被指数分析香蒲的反射光谱，在以后的估算模型研究中，可以使用其他光谱变换方法、光谱参数或引入其他的植被指数与砷含量作相关性分析，提高相关系数，进一步找出对砷敏感的波段，在建立模型方面可以通过随机森林、支持向量机等非参数模型，以期构建出估算效果更好的模型。