淮北矿区有机质含量反演

汤 超

(安徽理工大学空间信息与测绘工程学院，安徽 淮南 232000)

引言

精准农业作为我国目前的农业发展重要方向，是把农业和现代技术结合的一种现代化农业应用系统[1]。快速获取土壤理化性质成为实现精确农业的重要环节，因此光谱分析掌握土壤有机质含量对精准农业有重要意义。土壤有机质中含有各种碳类有机化合物，其中含有植物生长所必需的各种元素[2]，同时是地下微生物活动的生命之源，是把持土壤肥力和养分供应能力的重要因素[3]。随着工业化、城镇化的快速发展，矿山类型企业成为国民经济稳步增长的重要支点，伴随共生的重金属污染、生态环境恶化等问题，致使收集土壤有机质含量的变化进行分析显现得尤为重要。目前，高光谱遥感是较为成熟的现代技术用以快速监测土壤有机质含量的方式。

现行条件下，有2种获取高光谱数据的方式，光谱分析法和传统的实验室分析法。实验室法会挑选采集并添加标签，风干，筛选出杂物，最后通过化学方法如光度法、重铬酸钾容量法等方法计算有机质的含量。传统方式的野外地面实测具有建模精度高、速度快等优点，但是利用人工测量不可避免会出现实用区域小、周期长，且模型的延展性低等缺点[4,5]。高光谱遥感法，利用遥感技术中的高光谱获取范围广、周期短、实时性强，适用于大面积建模等优点，通过辨别地物之间细微的吸收和反射特性来实测地物的地球化学元素含量。从而得到光谱反射率和含量之间的曲线，通过分析地球化学元素含量的相关性，找到化学元素的敏感波段，建立线性或者非线性模型，以此反演化学元素的含量分布[6-10]。王婷等[6]通过对华南地区75组高光谱数据的选取，选择敏感的单波段和双组合波段，对比分析了线性回归和非线性的回归效果，得到最佳的精度模型，相对误差也更小。Lijuan Cui等[7]根据敏感波段用多元线性回归法建立有机质与反射率的模型，经过精度评价，其均方误差(RMSE)的值在0.38%～1.67%；全文诠释了有机质的建立和形成的长期变化，并建立模型[11-13]。目前有机质的研究，关注点在构建模型和光谱指标变换。因此，本文通过ASD便携式光谱仪对朱庄矿采集的样本进行密室探照，检测测定土壤光谱数据后用实验室分析法获取土壤具体有机质含量，构建安徽省淮北市朱庄矿区的土壤有机质高光谱估算模型，用建立的模型估测矿区的有机质含量，为今后的研究提供坚实的基础，进一步为淮北矿区采动裂隙评价和矿区土地复垦以及生态资环修复打下基础。

1 材料与方法

1.1 研究区概况

研究区朱庄煤矿位于安徽省淮北市中部，具体位于离淮北市区9km的安徽省北市杜集区矿山集镇境内。地处N33°16′～34°10′，E116°24′～117°03′。近年来，煤矿开采，矿区已多处产生裂隙。由于该地开采严重，生态系统的恢复能力较差，所以对该地进行有机质含量的反演有重要的实际意义。在朱庄矿将土壤有机质含量样品送往安徽理工大学分析测试中心，采用元素分析仪法测定土壤中含碳量，将测定的有机碳含量转化为有机质含量。

1.2 数据来源及预处理

1.2.1 土壤光谱反射率测定及处理

采用ASD便携式地物光谱仪对研究区土壤样本的光谱反射率进行了测量，有效的光谱波段为350～2500nm。在光源上选择与太阳光源相似的50W卤素灯，将土壤样品放入具有角度为10°探头的密室内测定以避免周围其它因素的干扰，在矿区的南北6个深度，取10个点位为样品点，用ASD viewSpecpro软件进行预处理，剔除异常曲线，在每个点位求平均值确定为每条光谱的反射率。

1.2.2 光谱数据分析

在MATLAB软件中采用Savitzky-Golay平滑方法，目的是在剔除异常曲线后将波形中噪声引起的细小突出平滑。如图1所示，经过处理后的光谱曲线分别在紫外可见光(350～1100nm)呈现上升趋势，但是在1300nm中有一个降低，之后又稳步上升。近红外(1100～2500nm)上在1500nm处有一个明显的凹陷，之后又缓步上升，但是上升趋势减缓。所有曲线的趋势大致相同，总体都呈快速上升，然后在1300nm和1600nm段有个呈现波谷式的下降，之后上升趋势减缓，推测是因为矿区粘土中所含的水分子或者羟基的吸收带。

图1 土壤反射率平滑后光谱曲线

1.2.3 朱庄矿土壤属性统计特征分析

土壤的理化性质包括土壤pH、有机质、有效磷、速效钾等各方面。而对土壤理化属性的描述统计方法包括变异系数、标准差、平均值、最小值、最大值等。平均值和均值是数理统计中广泛用于描述频率分布以及中心位置集中趋势的2种方法。变异系数(CV)用来反映数据集的离散程度。本文有机质含量变异系数为0.603(0.1≤CV≤1)，表现为中等变异性。变异系数越大，表明土壤理化性质中的空间离散程度越大，而变异系数越大，也有利于模型的构建，使模型更具普适性。

表1 土壤有机质特征分析

1.2.4 光谱变换

Al-Abbs[15]等研究表明，土壤有机质的含量与反射光谱呈负相关。在可见光与红外光的波段发现本实验数据满足以上关系[15-17]。根据以上关系，测定数据时可以根据有机质含量的起伏变化关系发现敏感波段，为建模打下基础[17]。在选择光谱变换指标时不应该只选择线性的，受到土壤和测量环境的影响，也应该考虑非线性的指标对反射率的变化。光谱的一阶微分和二阶微分能变化去除一些线性噪声和基线漂移等产生的影响，而进行吸光率变换(倒数的导数)也是处理方式之一，吸光率变换能解决一些非线性和光照条件的影响，增强光谱数据和土壤有机质之间的相关性。综上，选择反射率、一阶微分、二阶微分导数的对数变换作为光谱数据构建的指标进行处理，找到波形上剧烈起伏大且对土壤有机质含量敏感的波段。

表2 光谱指标描述

1.2.5 模型建立与验证

在筛选光谱指标上，利用皮尔逊(Pearson)相关分析方法进行筛选。光谱指标中的函数自变量选择反射率指标和经过变换后的反射率指标，因变量则选择通过实验室密室测得的土壤有机质含量，选择合适的模型进行回归反演。相关系数用于诠释事物或者变量之间相关性程度强弱关系，皮尔逊相关系数是最常见但不是唯一的相关系数。不同的事物或者变量之间对应着不同特征的统计指标。在光谱分析中，会分析许多变量之间的关系问题，如在土壤有机质反演中，会分析土壤有机质含量和反射率及其反射率变换的指标(曲线的为非线性)的相关关系，会分析土壤理化性质之间的相关关系，会分析土壤训练集与验证的得到的模型预测值与实测值之间的相关关系。会根据以上找到合适的相关性强的敏感波段，建立模型揭示土壤有机质和光谱的关系。而土壤预测值与实测值之间的关系一般条件下选择y=x模型检验最后的模型精度。

光谱分析中2个事物或者变量之间有正相关、负相关和无关等，本文分析高光谱有机质含量和光谱指标的关系，研究表明都呈现负相关。而2个变量的线性相关用相关系数r2或r表示。相关系数r公式：

(1)

2 相关性分析和结果

2.1 光谱指标与土壤有机质含量的相关性

相关性分析可用于建立土壤有机质含量和反射率光谱指标之间的相关性强弱关系，也是较为常用的进行光谱反演建模的分析手段。如图2-5所示，将土壤原始光谱数据和转换后的一阶微分、二阶微分、吸光率变换(倒数的对数)与土壤有机质含量进行Pearson相关分析。在图2中观察经过平滑后的反射率曲线与有机质含量的Pearson相关性变化趋势，正如大部分研究所示一样，有机质含量与平滑后的反射率曲线都是负数，呈现为负相关，且相关性程度随着可见光到红外光的波段呈降低趋势。在350～1761nm波段都表现出很好的相关性，在558nm处达到相关最大。一阶微分变换能够去除部分线性的噪声，增强局部位置的光谱响应差异。一阶微分变换数据与有机质含量在350～2500nm波段上可以明显地找到更多的吸收反射特征点，并在350～1760nm的波段上可以发现更多的拐点，在一阶微分变换的数据中，507nm、519nm、1399nm、1432nm的光谱波段表现十分显著。而在二阶微分中，二阶微分能够消除一些基线漂移和一些仪器背景的影响，二阶微分表现的光谱数据中，在350～1310nm段相关性起伏皆不大，1390nm、1399nm、1409nm、1446nm波段分别达到一个极大值或者极小值点，并且后面又趋于稳定。收集样本点的光谱信号时，由于土壤表面的散射和吸收作用引起吸收强度的变化造成一些非线性的光谱噪声。将反射率转化为吸光率是处理非线性的常用方法，将土壤反射率转化为吸光率不仅增强了光谱反射的差异性，还可以减少光照条件等其它随机因素的影响，更好地反映土壤的理化性质。吸光率和反射率的趋势相似，在极值点531nm段达到极大值点后，便不断下降。

图2 反射率与有机质含量的相关系数

图3 反射率一阶微分变换与有机质的相关系数

图4 反射率二阶微分与有机质的相关系数

由图2-5可知，土壤有机质反射率的光谱曲线的波形在可见光和红外光的波谱上除了羟基和水分子的影响都较为平直，明显有机质含量的吸收特征都不明显，所以需要将土壤光谱进行光谱指标变换。而变换之后的曲线，能扩大光谱特征中差异不大的部分，增大光谱的吸收特征，突出对有机质含量中较敏感的波段，方便识别各种光谱指标中敏感波段中的差异。将原始光谱反射率进行光谱一阶微分变换、二阶微分变换、吸光率(倒数的导数)转换等处理，从图2至图5可以看出,土壤光谱曲线在经过微分变换后，可见光与红外光波段上的光谱指标的细小差异被凸显出来。

图5 反射率吸光率变换(倒数的对数)与有机质的相关系数

2.2 土壤有机质含量的估算反演模型

偏最小二乘回归(PLSR)是一种主成分分析和相关分析相结合的多元统计回归方法。偏最小二乘回归的优点是当多个自变量相关时，更容易识别噪声和处理信息，且模型将包括所有自变量。本研究认为，PLSR模型更适用于土壤光谱估计分析，通过建立X(土壤光谱数据)与Y(土壤有机质含量)的最优模型，达到土壤光谱信息准确估算理化参量的目的。

在参与回归反演的PLSR模型，选择参与回归的变量个数非常关键。如果变量合格数太少，则模型的拟合程度不高；反之变量个数太多，则会出现过度拟合现象。所以选择交叉检验方法(K-CV)作为检验模型的精度和估算结果。在研究中，原始数据集被随机分成K个子样本，选择其中一个样本为验证数据集，剩下的K-1个样本用于训练集，最后总的重复验证K次。根据数据选择K的值为3，原始数据集被随机分成3个相等的组，其中任意2组为训练集，剩下的一组为验证机，重复3次上述操作。

2.3 精度评价

建模以后需要对建模进行精度检验，通常选择表现线性拟合程度的决定系数(R2)和表示观测值与真值偏差差异的均方根误差(RMSE)。平均相对误差(MRE)来预测模拟值与真值的精度。验证模型选用的指标是决定系数(R2)、平均相对误差(MRE)和均方根误差(RMSE)，计算公式分别如式(2)、(3)、(4):

(2)

(3)

(4)

3 讨论

为了消除噪声和基线漂移以及非线性噪声影响，将反射率进行4种光谱指标变换处理，分别为经过Savitzky-Golay平滑后的反射率曲线，一阶微分、二阶微分、吸光率(倒数的导数)变换。分别筛选出4种光谱指标中相关性较高也产生起伏的敏感波段。如表3所示。

表3 有机质含量的敏感波段及其相关系数

通过上述的相关性分析成果和交换检验方法，在SPSS软件中用PLSR建立模型，分别对不同的建模集和验证集进行分析，从而选择最佳的模型预测建模模拟预测值。土壤有机质PLSR的估算模型结果如表4所示。

表4 土壤有机质偏最小二乘法估算模型

由表4可知，在所有的反射率以及光谱变换形式下，反射率转换为吸光率具有较好的反演能力，其中吸光率无论是在建模集还是验证集R2,RMSE的值都有较好的准确性，R2值越大，RMSE值越小，说明吸光率变换有着很好的反映能力。在建模集的R2值达到0.45，在验证集达到最大0.657。在RMSE的值在吸光率变换时都很小，如在建模集为1.082，而在验证集达到0.932，都在很大程度上确保了吸光率变换的可靠性。吸光率模型：

(5)

4 结论

建模完成后对建立的模型检验可靠性，将验证集中有机质含量的实测数据与预测数值进行对比分析。模拟实测值与预测值的含量值在0.3～5.3g·kg-1，平均值是2.55g·kg-1，相关系数(R2)达到0.65，计算最后以反射率转化为吸光率的平均相对误差(MRE)为78.6%。上述数据和图像皆可说明土壤有机质用光谱吸光率来反演较为可靠。

图6 土壤有机质模型实测值与预测模拟值之间的比较

安徽省淮北市朱庄矿区开采多年，空间变异显著，该区域适合进行有机质含量反演。

光谱反射率与有机质含量呈现负相关，在紫外到可见光以及近红外都与有机质含量表现较好的相关性。其中与吸光率的变换形式最佳，可在该矿区推广使用。研究区的有机质反演模型以吸光率为自变量时，决定系数(R2)能达到0.45，均方根误差(RMSE)能达到1.08。当用吸光率模型预测有机质含量时，与实测值的决定系数(R2)达到0.65，平均相对误差(MRE)为78.6%。以上均说明预测模型精度高，稳定性好。

研究区为矿区，开采严重，出现地面水下跌、错动下沉、矿区裂隙等问题，研究有机质含量反演结果可为矿区生态修复、矿区复垦，以及矿区土地修复提供数据。

本研究对有机质含量反演的影响，在光谱变换、建模方法等方面存在缺陷，后面会考虑非线性的模拟方法，如决策时、人工BP神经网络、光谱变换指标选择小波变换等。