徐良骥，刘曙光，孟雪莹，韦 任

(1.深部煤矿采动响应与灾害防控国家重点实验，安徽 淮南 232001; 2.安徽理工大学 测绘学院，安徽 淮南 232001)

近年来矿区沉陷水域的研究主要集中在水体的氮、磷含量及分布特征[1-4]、富营养化[1,5]、水化学特征及影响[6-7]、微动物群落结构及影响因子[8-9]、水质[10]、水环境评价和综合利用分析[11-12]等;随着对重金属污染的重视，很多学者研究了海洋[13]、湖泊[14]、矿区附近河水[15-19]、沉陷水域[20-21]的重金属含量及分布特征、重金属污染评估、有害元素的健康风险评估[22]、重金属对微生物群落及生态环境的影响[9]。但关于沉陷水域中重金属含量反演模型的建立很少。本工作利用地物光谱仪采集沉陷水域的光谱反射率并对光谱数据进行微分变换，采用单波段回归法、逐步回归法和波段深度与偏最小二乘结合3种方法分别建立各重金属元素含量的预测模型，并对3种回归模型进行精度评定，确定各重金属含量的最佳预测模型。

1 研究区概况

潘一矿沉陷水域位于安徽省淮南市潘集区(地理坐标为北纬32.80°～32.83°，东经116.78°～116.82°)，该沉陷水域附近存在煤矸石堆、煤矸石充填复垦的试验田和动物养殖，该研究区包括两块沉陷水域，左侧沉陷水域目前用于水产养殖，右侧水域暂无工农业用途。研究区水域面积约为0.35 km2，共有43个采样点，采样点均匀分布，各采样点间距约为150 m(图1)。

图1 采样点分布示意

2 样品采集与分析

2.1 水体样品

研究水域内采集43个0～10 cm表层水体样品，每个样品取约2 L。采用原子吸收分光光度法测定水样中Cu,Zn,Pb,Cd和Cr五种重金属元素的含量，采用原子荧光光度法测定水样中As元素的含量，测定结果见表1。

2.2 光谱数据的采集

采用美国ASD公司生产的地物光谱仪FieldSpec 4在沉陷水域进行光谱采集，该仪器的光谱数据波段范围是350～2 500 nm，有两种采样间分别为1.4 nm(350～1 000 nm)、2 nm(1 000～2 500 nm);该实验的重采样的间隔为1 nm[23]。光谱数据为2017-07-27T10:00—12:30 在潘一沉陷水域进行采集，每个采样点采集20条光谱曲线，并将利用ViewSpecPro软件剔除掉异常曲线后的反射率的算数平均值作为样本的原始反射率光谱值[24]。

2.3 光谱数据预处

水体光谱数据采集过程中受随机误差影响，且重金属光谱响应信号微弱、原始光谱数据难以直接反应特征波段的特点，对原始光谱数据进行一阶微分变换、二阶微分变换和倒数对数变换，增强光谱响应特征波段，变换后的光谱曲线如图2所示。微分变换是常见的光谱数据处理方法。可在某种程度上消除、削弱光谱数据中的噪声，放大光谱信息，改善多重共线性，提高标定方程的性能[25]。计算公式为

(1)

(2)

表1 水样的重金属元素含量统计特征

Table 1 Statistical characteristics of heavy metal elements in water samples

编号元素含量/(mg·L-1)Cu2+Pb2+Zn2+AsCd2+Cr6+PB10.00500.00700.05400.00400.00060.0040PB20.02400.01700.11700.06700.00150.0050PB30.01200.01500.15700.02900.00120.0060PB40.00300.00400.14600.01200.00010.0040PB50.00800.01000.06900.04500.00010.0040PB60.00400.02200.01100.00400.00130.0100PB70.00100.01200.00600.03100.00120.0040PB80.00100.00400.05200.01600.00100.0030PB90.02000.01700.06200.03500.00100.0090PB100.00200.00100.04100.02000.00100.0010PB110.00600.00100.10500.00700.00100.0020PB120.00300.00700.03600.04100.00100.0010PB130.00100.02100.01700.00100.00100.0050PB140.00100.00700.00600.02400.00100.0010PB150.00500.00700.05400.00400.00060.0040PB160.02400.01700.11700.06700.00150.0050PB170.01200.01500.15700.02900.00120.0060PB180.00300.00400.14600.01200.00010.0040PB190.00800.01000.06900.04500.00010.0040PB200.00400.02200.01100.00400.00130.0100PB210.00100.01200.00600.03100.00120.0040PB220.00100.00400.05200.01600.00100.0030PB230.00790.01230.07510.03160.00060.0055PB240.00860.01050.08540.01930.00080.0049PB250.00600.00100.10500.00700.00100.0021PB260.00410.00870.09230.01850.00050.0050PB270.00740.01400.05930.01940.00090.0062PB280.01020.01430.07190.03410.00090.0058PB290.02000.01700.06200.03510.00100.0091PB300.00370.01070.05730.02710.00080.0049PB310.00210.01140.03440.02890.00100.0042PB320.00200.00100.04100.02000.00100.0010PB330.00820.00890.06000.02330.00100.0049PB340.00290.00700.03600.04110.00100.0010PB350.00100.02100.01700.00100.00100.0051PB360.00350.00910.08070.02090.00060.0046PB370.00520.00640.06030.01660.00090.0034PB380.00660.00670.08260.01760.00100.0033PB390.00960.01130.08510.02650.00110.0040PB400.00620.01160.04710.03000.00110.0032PB410.00640.00710.05830.02440.00100.0032PB420.00410.01380.05050.01610.00100.0041PB430.00100.00700.00600.02400.00100.0010

图2 沉陷水域光谱曲线

在对光谱研究时，将光谱原始的反射率进行倒数对数变换是常用的处理方法，该变换形式可以增强相似光谱间差异，减少因光照条件、环境变化等因素的影响，计算公式为

(3)

式中，λi+1,λi,λi-1为相邻的波长;R(λi+1),R(λi),R(λi-1)为相应波长的原始反射率;R′(λi),R″(λi),Rlgλi分别为波长λi对应的一阶微分、二阶微分、倒数对数光谱反射率值。

2.4 统计分析

采用决定系数(Coefficient of Determination,R2)、均方根误差(Root Mean Square Error，RMSE)对模型的可靠性以及它的预测能力进行评价[26-27]。R2越接近1，表示模型的预测值越接近于实测值，模拟精度越高。RMSE越接近0，说明拟合能力好。

3 结果与分析

3.1 重金属含量与水体光谱反射率的相关性

建立沉陷水域重金属含量反演模型的首要步骤是统计分析水体中重金属元素含量与水体光谱反射率数据的相关性。分析两类数值变量的相关性的常用指标为Pearson相关系数，计算公式为

(4)

式中，ri为光谱指标反射率与水体中重金属元素含量的相关系数;i为波段序号;cov(Xi,Y)为第i个波段水体光谱指标反射率和重金属含量间协方差;var(Xi)为第i个波段光谱指标反射率值的方差;var(Y)为水体样本重金属含量的方差[26]。

将6种重金属(Cu,Pb,Zn,As,Cd和Cr)含量分别与原始反射光谱(REF)、一阶微分光谱(FDR)、二阶微分光谱(SDR)、倒数对数光谱(lg(1/R))进行相关性分析，相关性分析结果如图3所示。由图3可知:

(1)Pb,As,Cd元素与原始光谱相关性较低，相关系数<0.3;且Cd元素与原始光谱在400～977 nm呈负相关。Cr相关性最大值在1 071 nm处，最大相关系数为0.43，为低相关。Cu,Zn元素含量与原始光谱相关性达到中度相关，分别在波段1 077,1 085 nm处出现最大相关系数为0.615 4,0.53。

(2)Cu元素、Cr元素与一阶微分光谱的大多数波段含量呈显著相关，相关系数波动范围为-0.768～0.874;Zn元素含量与部分波段达到中度相关，相关系数波动范围为-0.746～0.750;Pb元素、As元素含量与一阶微分光谱呈中度相关，相关系数波动范围为-0.444～0.566;Cd元素含量与一阶微分光谱呈低度相关，相关系数波动范围为-0.49～0.47。相较于原始光谱的相关系数，都有不同程度的提高。

(3)二阶微分光谱绝大多数波段与Cu元素、Zn元素、Cr元素含量达到高度相关，相关系数在-0.885～0.824;部分波段与Pb元素、As元素、Cd元素含量呈中度相关，相关系数在-0.648～0.764。

(4)倒数对数光谱与各金属元素相关系数较小，与原始光谱的相关性曲线接近且近似呈镜面对称分布。

3.2 基于单波段的水体重金属含量估算

选出29个水样进行建模，14个水样进行验证。选取400～1 200 nm波段范围内的光谱(REF,FDR,SDR, lg(1/R))与重金属含量相关系数最大的波段建立重金属含量反演模型，各重金属的最优拟合模型方程(表2)。

图3 光谱指标反射率与重金属元素含量相关性

表2 特征单波段与重金属含量最优拟合模型

注:拟合模型方程中，Y为水体重金属元素含量;x为水体光谱拟合波段的波长位置的反射率值(下同)。

相比于REF,FDR,lg(1/R)，SDR与各重金属含量的拟合模型方程效果较好。由表2可知，SDR与Cu,Zn,As和Cr元素含量的二次方程的R2均>0.5;该拟合方程拟合效果较好;SDR与重金属Cu含量的二次拟合模型精度最高，建模R2为0.796，REMS为0.000 28 μg/L;模型验证R2为0.823，REMS为9.66×10-6mg/L。Pb和Cd元素含量的建模R2均<0.5，验证R2也不高，表明该拟合方程拟合效果不好。

3.3 基于SMLR的水体重金属含量反演回归分析

将400～1 200 nm波段内的波峰、波谷和(在 0.05 检验水平上)与水体重金属含量成显著相关的光谱波段选为预测模型的备选波段。利用 SPSS 软件建立拟合模型，得到的拟合模型方程见表3。

由表3可知:同一金属元素的REF,FDR,SDR,lg(1/R)指标下，SDR与重金属(Cu,Pb,Zn,As和Cr)含量的拟合模型方程的建模R2均>0.7，表明SDR的特征波段建立的拟合方程能较好的拟合光谱数据与重金属含量之间的关系。Pb,Zn,As和Cd四种元素只有二阶微分光谱能与之建立较好的拟合方程。

表3中方差扩大因子VIF均<10，表明各多元方程的特征波段之间不存在多重共线性。

表3 水体重金属含量逐步回归模型

Table 3 SMLR model for heavy metal content in water

光谱指标光谱波段/nm拟合模型方程建模R2RMSE/(mg·L-1)R2验证RMSE/(mg·L-1)共线性统计容差VIFCuREFFDRSDRlg(1/R)116110771158672719785100511671076Y=-0.070x1161+0.091x1077-0.0030.5058.97×10-60.5431.44×10-5Y=-11.653x1158+0.0010.1201.64×10-50.1372.81×10-5Y=107.217x719-35.718x672-45.799x785-107.217x1005+0.0140.7361.22×10-40.7701.36×10-4Y=0.011x1167-0.035x1076+0.0230.5119.08×10-60.5191.57×10-50.1536.5550.1536.5551.0001.0000.4822.0740.4332.3100.6061.6510.8341.1990.1516.6330.1516.633PbFDRSDR9985658641045Y=7.32x998+0.0050.1732.16×10-50.2103.70×10-5Y=-146.102x565-198.881x864-277.202x1045+0.0010.7456.6×10-60.7721.07×10-51.0001.0000.3313.0240.2264.4250.4502.224ZnREFFDRSDRlg(1/R)10851153687845101311654794836831074Y=0.31x1085+0.0180.2629.46×10-40.2811.72×10-3Y=0.162x1153+0.0560.1910.0013640.3932.42×10-3Y=-5.901x479-501.493x483+116.69x683+433.8486x687-1270.793x845+384.686x1013+636.293x1165+0.1020.7752.89×10-40.7745.37×10-4Y=-0.101x1074+0.1550.2279.91×10-40.2201.86×10-31.0001.0001.0001.0000.2434.1180.3442.9040.4832.0680.5861.7070.2813.5580.2743.6450.2833.5361.0001.000AsFDRSDR418411704Y=27.977x418+0.0270.2861.26×10-40.3202.26×10-4Y=79.1021x411+751.033x704+0.0160.7095.10×10-50.7249.12×10-51.0001.0001.0001.0001.0001.000CdFDRSDR42842298710511138Y=-0.783x428+0.0010.2901.60×10-70.2411.97×10-7Y=0.735x422-2.424x987-8.974x1051-5.251x1138+0.0010.5654.05×10-80.5646.70×10-81.0001.0000.5291.8900.5571.7950.9551.0480.5921.689CrFDRSDR604657566Y=5.433x604+8.025x657+0.0030.7669.91×10-70.7671.68×10-6Y=-28.710x566+0.0030.7811.07×10-60.9021.67×10-60.5961.6770.5961.6771.0001.000

3.4 波段深度与PLSR结合的水体重金属含量反演回归分析

连续统去除后计算如下光谱吸收特征:BD(波段深度)、BDR(波段深度比)、NBD(归一化面积波段深度)、BNA(归一化面积波段指数)，表达式为

BD=1-R

(5)

BDR=BD/BDmax

(6)

NBDI=(BD-BDmax)/(BD+BDmax)

(7)

BNA=BD/BDarea

(8)

R(连续统去除反射率)是光谱反射曲线上的波长反射率与对应连续统线上数值的比值。选取光谱的波段深度作为基础数据，计算出该光谱的特征指数，通过偏最小二乘回归方法构建与水体重金属含量的预测模型，模型精度见表4。

由表4可知:重金属Cu,Cr元素与BD参数结合PLSR建立的模型的决定系数均大于0.5，表明重金属Cu,Cr元素含量与该模型预测值具有较好的的线性关系。其中，BNA结合PLSR的模型的预测值与Cu,Cr元素含量的线性关系最好，决定系数分别为0.783,0.758，均方根误差为5×10-7,1.7×10-6mg/L。NBDI,BNA结合PLSR能建立Zn元素含量反演效果较好的模型;BDR,NBDI,BNA能建立As,Cd元素反演效果较好的模型;NBDI能建立Pb元素含量反演效果较好的模型。

4 结 论

(1)通过对沉陷水域的光谱数据进行不同形式的变换可知，一阶微分变换、二阶微分变换、倒数对数变换对光谱信息都进行不同程度的放大。① 从相关系数来看，一阶微分变换、二阶微分变换效果较为明显;微分变换后的光谱波段与重金属含量达到高相关度;Cu元素与一阶微分光谱的相关系数最大，为0.874;Pb元素、Zn元素、As元素、Cd元素和Cr元素与二阶微分光谱的相关系数最大，相关系数分别为-0.648,0.824,0.764,0.636,-0.885。② 根据单波段回归和多元逐步回归建模的决定系数来看，二阶微分光谱对6种重金属(Cu,Pb,Zn,As,Cd和Cr)元素含量的拟合效果效果最好;即二阶微分变换的预处理方法比较适合。

表4 波段深度结合偏最小二乘反演重金属结果

Table 4 BDA combined with PLSR inversion of heavy metal content model accuracy

指标BDBDRNBDIBNACuR20.6150.6920.7680.783RMSE/(mg·L-1)4.30×10-61.56×10-51.45×10-45.00×10-7PbR20.1560.3460.5090.465RMSE/(mg·L-1)5.80×10-64.50×10-64.56×10-57.00×10-7ZnR20.4860.2130.5580.643RMSE/(mg·L-1)3.70×10-64.30×10-66.90×10-65.00×10-7AsR20.2330.5590.6780.511RMSE/(mg·L-1)1.12×10-53.64×10-54.50×10-61.40×10-6CdR20.1590.5100.5120.499RMSE/(mg·L-1)1.55×10-54.53×10-51.85×10-57.00×10-7CrR20.5980.6100.7200.758RMSE/(mg·L-1)2.56×10-53.69×10-51.23×10-51.70×10-6

(2)根据单波段回归，多元逐步回归和波段深度与PLSR结合3种模型的R2和RMSE大小可知，二元微分的单波段回归对Cu,Cr元素含量的拟合效果最好;Cu元素建模R2为0.796，RMSE为2.8×10-7mg/L，模型验证的R2为0.823，RMSE为9.66×10-6mg/L;Cr元素建模R2为0.794，RMSE为5.63×10-7mg/L，模型验证的R2为0.806，RMSE为1.65×10-6mg/L;二元微分的多波段逐步回归对Pb,Zn,As和Cd四种元素的拟合效果最好，建模R2分别为0.745,0.775,0.709,0.565，RMSE分别为6.6×10-6,2.89×10-4,5.1×10-5,4.05×10-8mg/L;模型验证R2分别为0.772,0.774,0.724,0.564，RMSE分别为1.07×10-5，0.537，9.12×10-5，6.7×10-8mg/L。