基于高光谱的土壤碳酸钙含量估算模型研究

吴倩，姜琦刚，史鹏飞，张莉莉

(吉林大学地球探测科学与技术学院，长春 130026)

0 引言

土壤是自然界不可或缺的资源[1]，它为作物的生长提供了一定的养分和水分。碳酸钙(CaCO3)是土壤中碳酸盐的主要成分,对土壤物理、化学性质、碳的沉积、土壤酸碱度、土壤胶体性状等的研究都有重要意义[2-3]。近年来，随着遥感科学的发展，高光谱遥感应用的深入，人们越来越体会到定量遥感的必要性，基于高光谱的土壤理化性质快速、定量、准确的模型构建成为农业遥感研究领域研究的重点[4]。土壤的光谱反射特性是土壤基本特性之一，与土壤中的各种理化性质有着十分密切的关系，土壤中各种物质特殊的光谱特征以及土壤组成物质的多样性都会在光谱曲线上体现，碳酸钙的含量对土壤光谱反射特性具有一定的影响，Gaffey[5]研究发现碳酸盐矿物由于碳酸根的振动过程，在可见-近红外(350～2 550 nm )光谱区域有7个较强的吸收特征； 张婷婷[6]研究认为土壤中碳酸盐含量对土壤的光谱特征有一定的影响，在2 230～2 270 nm之间有明显的吸收特征； Stenberg[7]研究认为土壤碳酸盐在2 335 nm处有较强的吸收特征； 以往研究表明碳酸钙的可见光、近红外(VNIR，400～1 200 nm)、短波红外(SWIR，1 200～2 500 nm)光谱可作为土壤碳酸钙含量估算的理论依据，是一种物理性的无损伤和低成本的研究方法，可以用来区分土壤类型和评估大量土壤特性。为了更好地理解电磁波与地表特征之间的相互作用，国内外学者将经典的数学物理理论与遥感实践相结合，建立了近百种定量遥感模型，这些模型大体可分为统计模型、物理模型、半经验模型3种。利用高光谱数据估测土壤理化性质的模型有很多，近年来的研究中大多运用统计模型，建立遥感数据与地面观测数据之间的线性回归方程，统计模型简便易行，一般仅包含3～6个参数。以往研究中多运用偏最小二乘回归、多元线性同步回归、支持向量机等方法，且研究内容大多集中在对土壤有机质含量、土壤水分、土壤游离铁含量的估算建模上。本次研究选取随机森林回归算法(random forest regression，RFR)对土壤碳酸钙含量进行估算，随机森林回归是一种集成学习模型，其优化了组合机器学习和CART决策树算法，并且对异常值和有噪声的数据具有一定的容忍度[8]，目前更多的应用于数据分类和回归分析[9]。

研究选取的土壤类型为黄绵土，黄绵土质地疏松，呈浅棕色，剖面发育不明显。广泛分布于甘肃、陕西、山西等水土流失较严重地区，是陕西省分布面积最大的土壤类型，其肥力低下，主要分布在坡度大、植被较少地区，是造成土壤养分有效性低、pH高，影响土壤结构及加速土壤侵蚀的重要因素之一，因此基于高光谱的土壤碳酸钙含量估算可以准确快速地获取土壤属性信息以及为水土流失地区提供治理等的理论参考，对于评价成土因素以及对土壤性质的研究具有重要意义。

1 材料与方法

1.1 研究区概况

研究区主要位于陕西省北部的黄土高原地区。107°15′41″～111°14′40″E，34°12′21″～39°29′17″N，黄土高原属于干旱区与半干旱区； 在气候分区上，该地区横跨温带和暖温带，年平均降水量约在300～700 mm之间，长期的水土流失导致生态环境较脆弱。陕西省土壤类型众多，本次研究的土壤类型为黄土高原上分布面积最大的黄绵土，研究共采集了78个土样，其中28个位于陕西省最北部的榆林市，41个位于延安市，9个位于渭南市。图1 为研究区地理位置及采样点分布。

图1 研究区地理位置及采样点分布Fig.1 Location of the study area and the distribution of sampling sites

1.2 样本采集与分析

根据陕西省黄绵土空间分布特征，于2016年7—8月期间在榆林市、延安市、渭南市共采集具有典型代表性的78个黄绵土土壤样品，对土壤样本进行风干、碾磨及过1 mm孔筛等操作后，采用四分法，将每个样品分为2份，一份利用气量法进行实验室土壤碳酸钙含量测定，另一份用作土壤反射光谱测定。

1.3 光谱测定

土壤VNIR/SWIR反射光谱是在暗室中控制辐照度条件测量的，测定采用了美国 Spectra Vista 公司生产的SVC HR-1024i便携式光谱仪进行，其波长范围为350～2 500 nm。在0.6～1 m之间放置1盏50 W白炽灯，设置60°方位角和15°天顶角，这种配置确保了由于阴影造成的最小扰动，使测量结果不那么依赖于粗糙度。用8°透镜提供直径为0.13 m的光斑。测量中利用白板进行校正，获得了3个连续的和50个单独反射测量系列，噪声水平低于0.001，光源的照射方向与垂直照射之间的夹角为30°，测量土样在4个方向上的光谱曲线。

1.4 数据预处理

为了提高光谱数据的质量，研究利用ENVI5.3对建模前的数据进行光谱数据库的建立以及10 nm间隔的重采样，重采样处理后每个土壤样本的光谱数据曲线包含216个波段。重采样后的光谱曲线依然存在一定量噪声，为了提高预测模型的精度，采用Savitzky-Golay 卷积平滑法(S-G)对数据进行进一步的处理[10]，平滑去噪后，对土壤碳酸钙光谱数据进行数学变换，得到原始光谱(R)的一阶微分(First order differential，R″)数据，二阶微分(Second order differential，R″)数据，连续统去除(Continuum Removal，CR)数据。

1.5 敏感波段选取

采用相关分析法(analysis of correlation，CA)与连续投影算法(successive projections algorithm，SPA)分别进行土壤碳酸钙含量与光谱反射率的相关性分析以及敏感波段的筛选，并依据建模结果对两种方法进行评价。

CA是指对两个或多个具备相关性的变量元素进行分析，从而衡量两个变量因素的相关密切程度[10]。 需要将4种光谱处理方法与所测样本的碳酸钙含量的真实值逐个波长进行相关性分析，并进行p=0.01显著性检验，根据相关性的高低来确定敏感波段，相关性越高，波段越敏感。相关分析采用SPSS中的相关性函数进行，通过显著性检验的波段作为接下来的建模的敏感波段。相关系数公式如下：

(1)

SPA在多自变量校正领域是一种采用前向选取的特征波段变量选择的方法[12]。能够有效地减少信息重叠，使变量之间的共线性达到最小，大大减少了建模变量的个数，有效提高建模的效率。实现用较少的信息量代表多数样本的光谱信息的效果[10]，SPA以均方根误差(ROOT mean square error，RMSE)为评价指标，将RMSE最小值下的波长个数确定为敏感波段数[13-14]，连续投影算法在Matlab2018a中实现。

1.6 模型建立与精度验证方法

随机森林(RFR)是一种集成算法(Ensemble Learning)，这种算法优化了组合机器学习和CART决策树算法[15]。随机森林有很多的优点，首先它能在决策类别时，去评估变量的重要性，其次它对于数据集较多的输入变量的处理精度较好，随机森林的两个应用主要是分类与回归。近几年随机森林在回归分析方面的应用也越来越广泛。该算法的估算精度比较高，树形结构自由生长，能够有效避免其他估算模型中的过度拟合现象，是对传统的决策树模型的一种优化，对数据量较大的数据集比较适用[15]。随机森林回归模型的构建运用R.3.4.1中的Random Forest包进行。

研究选取土壤碳酸钙含量实测值与估测值的均方根误差(RMSE)，决定系数(R2)以及相对分析误差(residual predictive deviation，RPD)3个参数对模型进行精度验证，均方根误差与相对分析误差越小，决定系数越大，模型的精度越高。

2 结果与分析

2.1 土壤碳酸钙描述性统计

首先对全部78个土壤碳酸钙样本进行异常值分析，将样本中碳酸钙含量不在±3倍标准差中的极度异常的4个样本剔除，从剩余的74个样本中选择出52个样本用于建模，22个样本用于验证，表1为土壤样本碳酸钙含量的统计结果，样本总体碳酸钙含量范围为11.109～175.077 g·kg-1，其中建模集中的52个碳酸钙含量范围为11.109～173.509 g·kg-1，验证样本碳酸钙含量范围为11.236～175.077 g·kg-1，样本的选取基本上涵盖了总体样本中的数据范围。建模集与验证集的偏度为0.46和0.25，均在0值附近，基本服从正态分布，建模集、验证集与总体样本的偏度差异不大，样本适合用于建立模型。

表1 土壤样本碳酸钙含量的统计结果Tab.1 Statistics of calcium carbonate in soil samples

2.2 土壤碳酸钙含量与光谱相关性分析

图2(a)为全部碳酸钙样本波长与反射率之间的原始光谱图，由图中可以看出，不同碳酸钙含量的土壤反射率曲线总体趋势基本相同，深度有所不同，吸收特征基本相同，在1 400 nm和1 900 nm和2 200 nm处有明显的吸收特征，这与前人分析结果相似，通常认为是由OH，AL-OH等引起的[16-17]，除此之外，近红外波段的反射率总体大于可见光部分的反射率； 如图2(b)所示，随着碳酸钙含量的增加，光谱反射率也随之增加，两者呈现正相关态势，且近红外波段的增长速度明显高于可见光部分的增长速度，反射率曲线在于2 350 nm处有一吸收谷，这主要是由于碳酸根的影响作用引起的[16]。

图3为基于相关分析算法的敏感选波段筛选图，由图3可以看出，原始光谱曲线在全波段的相关性较差，相关系数最高值在2 340 nm左右，相关系数仅为0.37，曲线整体呈现正相关的态势。经过一阶微分、二阶微分、连续统去除后，光谱曲线与土壤碳酸钙含量的相关性得到了显著的提高，从而可以进一步提高模型的预测能力。其中一阶微分相关系数最大值在1 610 nm处达到了-0.66，二阶微分相关系数在2 230 nm处达到了0.56，微分运算虽然可以较好地消除背景信号，但高阶微分在应用中往往表现出不稳定性，在应用中存在一定的问题，以往研究发现，在光谱分析中随着阶数的增加，相关系数呈现先增加后减少的趋势，并在二阶处达到最大值[18]，因此，研究采用一阶与二阶微分。连续统去除相关系数在2 150 nm处为到0.53，均达到中度相关水平，4种数学处理方法分别有54，82，54，44个波段入选建模波段。

图4为基于连续投影算法的敏感波段筛选图，经过连续投影算法降维处理过后筛选的敏感波段在图中用方块圈出，经S-G平滑后的原始光谱数据，RMSE最小值为17.966 8，共350 nm，440 nm，480 nm，620 nm，770 nm，840 nm，920 nm，970 nm，1 670 nm和1 880 nm等17个波段入选； 经一阶微分后的光谱数据，RMSE最小值为26.194 5，共420 nm，470 nm，620 nm，1 510 nm和1 530 nm等9个波段入选； 经二阶微分后的光谱数据，RMSE最小值为21.272 4，共360 nm，420 nm，830 nm，880 nm和1 020 nm等16个波段入选； 经连续统去除后的光谱数据，RMSE最小值为22.853 3，共410～430 nm，540 nm和570 nm等15个波段入选。与相关分析法相比大大减少了模型的输入变量，可有效提高模型的计算速度。

2.3 建模结果分析

表2为土壤碳酸钙含量敏感波段RFR模型，研究表明，经数学变换的光谱为变量建立的模型精度明显高于原始光谱，数学变换有效提高了土壤碳酸钙含量与光谱反射率之间的相关性，验证集R2均大于0.58，其中二阶微分的估算精度最高，验证集R2为0.82，RPD为2.37，表明模型较为适用； 从敏感波段选取方法看，除基于CR的建模集精度相关分析法较高于连续投影算法外，其余数学变换方法无论建模集验证集基于连续投影算法建模精度均高于相关分析法，由此看出，连续投影算法不仅变量少效率高，并且在波段筛选方面较为适用； 就模型本身而言，随机森林回归整体的预测精度较好，验证了利用随机森林回归估算土壤碳酸钙含量的可行性。就模型的稳定性而言，建模集整体精度均高于验证集，模型缺乏一定的稳定性，可能原因为： ①土壤碳酸钙样品普遍存在标准差较大、变异性较大的现象，碳酸钙的含量较为分散，在对碳酸钙含量进行估算研究中数据普遍存在以上问题； ②获取敏感波段的方式不同。例如，洪长乔等[17]除了利用相关分析法选取的波段，还添加了前人研究结果较好的几个数据集作为自己研究的建模集之一，并分别建模进行比较； Gomez等[19]在研究拉比因河谷土壤碳酸钙含量时，利用可变重要性(VIP分数)和B系数的值来确定光谱带，由此可看出敏感波段的选取方式是多样的，选择方式的不同对模型精度与稳定性会产生一定的影响。图5为基于SPA的4种数学变换RFR模型实测值与预测值拟合图，可以明显看出二阶微分的模型拟合效果较好。

表2 土壤碳酸钙含量的敏感波段的RFR模型Tab.2 RSR models for soil calcium carbonate content based on sensitive bands

3 结论

研究采用随机森林回归的方法建立了陕西省北部黄绵土高光谱碳酸钙估算模型，旨在探讨不同光谱处理方法与不同敏感波段筛选方法模型构建的优劣以及利用随机森林回归估算碳酸钙含量的可行性，主要取得了以下3点结论：

1)碳酸钙较明显的吸收特征位于1 610 nm，2 220 nm，2 340 nm，2 270 nm等处，且土壤碳酸钙含量与全波段光谱呈现较明显的正相关关系。

2)基于连续投影算法的模型构建精度除建模集CR变换外精度均高于相关分析法，表明连续投影算法在波段筛选中较为适用。

3)基于连续投影算法与二阶微分变换的随机森林回归建模精度最高，建模集R2为0.89，验证集R2值为0.82，且RPD值大于2，模型可信度较高，验证了利用随机森林回归估算土壤碳酸钙含量的可行性。