东北黑土区主要分布在吉林省和黑龙江省, 少量分布在辽宁省北部和内蒙古东部, 如图1所示。 根据以往研究, 黑土区土壤肥沃, SOM含量丰富, 且自南向北递增, 极少部分SOM含量少于2%。 采样点位于松辽平原北部的黑土带内, 主要分布在黑龙江省中部, 海伦市(44° 04'—44° 45'N, 125° 26'—126° 31'E)和吉林省中部, 德惠市(47° 00'—47° 41'N, 126° 25'—127° 06'E)。

图1

Fig.1

	Figure Option ViewDownloadNew Window
	图1 研究区及采样点分布Fig.1 Study area and sampling regions

土壤样本采集于2017年5月4日至9日, 正值春耕时期, 土壤表面冰雪完全消融, 且无植被覆盖。 共有74个采样点, 每个采样点在5 m×5 m范围采集5份土壤, 采集深度为0~15 cm, 所有样本均被独立装入自封袋并标记样本号, 采样点中心地理坐标由便携式GPS记录。

在实验室条件下将每个采样点的5份土壤混合均匀, 过2 mm筛剔除碎石子残留秸秆等杂质。 每份样本共取出100 g土壤分别用于测定各土壤成分含量, 其中, 利用重力法测定SMC含量, 重铬酸钾加热法测定SOM含量, 点位法测定土壤pH值, 等离子光谱仪测定Fe含量。 土壤样本成分含量统计如表1所示。 其中SOM含量介于2.0%~6.3%, 平均值约为3.5%, 整体含量较高; 土壤样本SMC近似于野外土壤水分含量, 值域分布区间较大, 介于3.8%~24.4%; Fe含量介于2.3%~3.4%之间, 且标准偏差(SD)与变异系数(CV)均较小(SD=2.12, CV=7.81), 总体含量差异不明显; pH值介于4.7~8.1之间, 均值为5.7, 标准偏差较小(SD=0.77), 总体偏离平均值较小。

表1

Table 1

表1(Table 1)

表1 研究区土壤样本成分值统计 Table 1 Statistical values of soil component from study area Contents Max Min Mean SD CV/% SOM/(g· kg-1) 63.36 19.87 34.86 12.65 36.29 SMC/(g· kg-1) 244.61 38.02 105.25 46.82 44.54 Fe/(g· kg-1) 33.63 23.02 27.16 2.12 7.81 pH 8.09 4.73 5.74 0.77 13.41

表1 研究区土壤样本成分值统计 Table 1 Statistical values of soil component from study area

室内光谱测量在黑暗的实验室进行, 将过筛土壤样本置于深度1 cm的铝制测试盒, 刮平表面, 用90 W的卤素灯照向表面。 采用ASD FieldSpec3光谱仪测定光谱反射率, 每个样本在不同位置共计测量20条光谱曲线。 利用ViewSpecPro软件计算平均光谱作为该样本的原始反射光谱曲线。

光谱数据首先以原始光谱反射率为输入矩阵, 进行标准正态变量变换。 逐一将光谱反射率数据X=(X₁, X₂, …, X_p)在p波长处转换为Z=(Z₁, Z₂, …, Z_p), 其中

Zi=(Xi-m)/s(1)

式(1)中, m和s分别表示X₁, X₂, …, X_p的均值和标准差^[13]。 以变换后的数据Z_n×p(n为样本数)作为输入矩阵进行多尺度连续小波变换。 光谱反射率的连续小波变换可以表示为小波基函数Ψ (λ )经过比例缩放和平移所产生的一系列系数, 用式(2)表示为^[12]

Wf(a, b)=< f, Ψa, b> =∫-∞+∞f(λ)1aΨλ-badλ(2)

式(2)中, W_f(a, b)为小波系数, a和b均为实数, a为比例因子(小波基的宽度), b为平移因子(小波基的位置)。 小波系数W_f(a, b)为小波基与光谱反射率在不同尺度下的关系。

参与对比的数据处理方法包括一阶导数(first derivative, FD), 二阶导数(second derivative, SD), 连续同去除(continuum removal, CR), 标准正态变量变换(SNV), 多元散射矫正(MSC), 高斯滤波法(Gaussian filter, GS), 以及对数变换[logarithm transform, log(1/R)]。 其中FD和SD变换, 首先对离散光谱反射率进行具有11个平滑点的二阶多项式拟合, CR通过调用Tripack R包中的Delaunay三角测量的凸包实现光谱拟合; 其余处理方法均利用Python实现, 且均以湿土壤样本的原始光谱反射率为输入数据。 为减少反射率边缘的噪声, 原始光谱反射率仅保留400~2 500 nm范围的值。 为了可视化不同处理方法的差异, 所有的信号变换都绘制在图2中, 其中SNV-CWT处理方法绘制的是10尺度(2¹⁰)的分解系数。

图2

Fig.2

	Figure Option ViewDownloadNew Window
	图2 经过不同处理方法的光谱曲线Fig.2 Reflectance spectraunder different transformation methods

皮尔森相关系数(Pearson correlation coefficients, PCC)常用来衡量特征与目标变量间的线性相关程度, 定义为两个变量之间的协方差与标准差的商。 基于模型的系数(model based coefficients, MBC)是一种基于模型的相关关系指标, 与传统的线性方法相比, 可以在建模之前判断特征变量与目标变量间基于模型的非线性关系。 灰色关联分析(grey relation analysis, GRA)是一种多因素统计分析方法, 通过计算参考数列和比较数列变化曲线几何形状相似程度判断因素间的关联程度, 由此方法计算得出的灰色关联度(GRD)可作为判定特征数据与相关因素序列的相关性指标。 研究中采用三种方法对SNV-CWT处理前后的光谱数据进行特征分析, 判断不同土壤属性的响应波段范围。

随机森林(random forest, RF)回归是一种基于决策树的机器学习算法, 采用bootstrap抽样技术抽取训练集, 并建立每个训练集的CART树, 构成森林, 具有强大的预测能力, 已被广泛应用于不同科学领域。 该方法学习速度快, 可以处理大量的输入变量, 即使在特征变量数远大于样本数据时也可以很好地避免过拟合, 且模型抗粗差能力较强。

采用决定系数(coefficient of determination, R²)、 均方差(mean squared error, Mse)和相对分析误差, 即标准差与预测均方根误差比(residual predictive deviation, RPD)这3个指标检测与评价模型的稳定性和预测能力。 其中R²越接近1说明模型拟合程度越好; 当RPD< 1.4时, 模型质量较差, 无法对样本进行预测; 当1.4< RPD< 2时, 模型质量较好, 可以用来进行目标变量的粗略估算; 当2< RPD< 2.5时, 模型质量很好, 可以用于目标变量的定量预测; RPD> 2.5时, 模型具有极好的预测能力^[14]。 检验与校正模型的统计指标(statistical indicators, SI)计算公式如式(3)—式(5)

Mse=∑i=1n[F(i)-A(i)]2n(3)

R2=1-∑i=1N(Fi-Ai)2∑i=1N(Fi-A̅)2(4)

RPD=SDMse

其中, SD=∑i=1n(F(i)-A̅)2n-1(5)

其中, SD为标准差, F(i)为预测值, A(i)为测量值, A̅为测量值的均值, n为样本数。

原始光谱反射率曲线在1 400, 1 900及2 200 nm附近均具有明显的拐点。 湿土壤条件下的光谱反射率较干土相比, 在1 400和1 900 nm附近具有较大程度的波动。 如图3所示, 两条曲线的有机质含量相近, 分别为6.0%和6.1%, 但含水量差异较大, 分别为15.8%和5.4%。 干土情况下两样本的光谱反射率曲线形状相似, 含水量不同时, 湿土样本两条光谱反射率形状具有较大差异, 说明水分对反射率具有显著影响。 而Fe、 pH与反射率的最大PCC分别为0.18与0.25, 与SOM的最大值(0.63)相比对反射率的影响较低。

图3

Fig.3

	Figure Option ViewDownloadNew Window
	图3 干湿土壤样本原始光谱反射率示意图Fig.3 Raw reflectance spectra of dry and wet samples

以建模精度为标准验证数据处理方法的有效性, 将置信度小于0.05的多维变量作为RF模型输入特征, 8种处理方法的预测精度如表2所示。 除log(1/R)变换外, SOM的建模精度均有一定的提高, 且经过SNV-CWT变换的模型精度可以达到R²=0.90, Mse=0.15%, PRD=3.21。 SMC在原始反射率情况下, 没有符合显著性水平的波段, 经过不同方法处理后仅GS与log(1/R)方法下的建模精度RPD< 2, 其余6种处理变换均有一定的提高, 且经过SNV-CWT变换精度为R²=0.93, Mse=1.37%, PRD=4.71。 常用方法中SD, CR, SNV, MSC和GS可以粗略建立Fe的预测模型, 然而建模结果精度均不理想, 经过SNV-CWT处理后建模精度可以得到最大的提升, 但RPD仍小于1.5。 8种处理方法较原始反射率相比, pH建模精度均有不同程度的提高, 且经过CR, SNV和SNV-CWT处理后, 模型精度RPD均可大于1.5。 由表2 可以看出, 不同数据处理方法对每种成分的影响效果不一, 采用SNV-CWT处理, 可同时获得四种土壤成分大幅度提升的建模精度。 与Stenberg等研究结果相比, SOM和SMC建模精度一般在 RCV2=0.50到0.90之间^[4], 经过SNV-CWT处理后, 在以多维特征潍输入变量的情况下SOM与SMC的建模精度已达到可获得的较好结果(R²=0.90和0.93)。

表2

Table 2

表2(Table 2)

表2 不同数据处理方法下的SOM, SMC, Fe和pH验证模型精度 Table 2 Statistical indicators ofvalidation models of SOM, SMC, Fe and pHunder different transformation methods Contents SI Raw FD SD CR SNV MSC GS Log(1/R) SNV-CWT R2 0.43 0.83 0.70 0.73 0.84 0.68 0.74 0.43 0.90 SOM Mse 0.83 0.25 0.44 0.40 0.23 0.47 0.38 0.83 0.15 RPD 1.46 2.51 1.85 2.01 2.71 1.76 1.97 1.37 3.21 R2 - 0.83 0.87 0.77 0.79 0.74 0.36 0.02 0.93 SMC Mse - 3.09 2.32 4.28 3.83 4.831 1.72 18.68 1.37 RPD - 2.79 3.49 2.45 2.34 2.34 1.26 1.13 4.71 R2 - - 0.20 0.17 0.14 0.09 0.27 - 0.48 Fe Mse - - 0.03 0.03 0.03 0.03 0.03 - 0.02 RPD - - 1.12 1.10 1.11 1.05 1.18 - 1.39 R2 - 0.44 0.36 0.64 0.61 0.28 0.39 0.13 0.62 pH Mse - 0.27 0.32 0.18 0.19 0.36 0.30 0.88 0.19 RPD - 1.38 1.27 1.67 1.60 1.18 1.28 0.95 1.63

表2 不同数据处理方法下的SOM, SMC, Fe和pH验证模型精度 Table 2 Statistical indicators ofvalidation models of SOM, SMC, Fe and pHunder different transformation methods

研究中以Mexh为小波基函数, 进行10个尺度(2¹, 2², 2³, …, 2¹⁰)的SNV-CWT分解, 将一维光谱反射率分解为10个尺度的小波系数。 以原始光谱数据集74×2 101维矩阵为输入, 经SNV-CWT处理后转换为74×10×2 101维矩阵。 针对四种土壤成分, 选择每个尺度置信度小于0.05的小波系数作为输入变量, 建立RF回归模型, 以验证模型获得的决策系数R²为标准, 确定各土壤成分最佳分解尺度。 如图4所示, SMC的模型精度在不同尺度上整体优于SOM, 且SOM和SMC的建模精度较Fe和pH更好。 这是由于SMC和SOM对光谱反射率影响起主导作用, 对Fe与pH验证模型信息具有掩盖作用, 且Fe的含量变化区间小。 SOM与SMC的最佳分解尺度分别为7和8; Fe与pH验证模型的最佳分解尺度为1和10。 每种土壤成分对应于不同的分解尺度有明显的增强信息, 相比于单一尺度分解, 多尺度SNV-CWT分解对多属性探测具有更大优势。

图4

Fig.4

	Figure Option ViewDownloadNew Window
	图4 SNV-CWT多尺度情况下的不同目标变量的模型决策系数Fig.4 Coefficients of determination (R2) corresponding to different soil properties on multi SNV-CWT scales

分别计算SOM, SMC, Fe和pH与原始反射率以及在最佳分解尺度系数的PCC, MBC和GRD, 其中PCC为线性相关系数的平方, MBC为单变量RF模型决策系数, GRD计算时取分辨系数ρ =0.5。 从图5中可以看出, 经过SNV-CWT变换后, 三种方法计算得出的相关系数均较原始反射率有所提高, 且SOM与SMC相关系数均高于Fe与pH模型。 这也充分解释了多维特征建模时Fe与pH的模型精度在不同尺度上均不及SOM与SMC。 与土壤成分对应的相关波段通过分析可分为两类: 一类是三种方法计算的相关系数曲线在相同的波长范围均有较明显的峰值, 这种波段可以认为是该成分的主要特征波段; 另一类是一种或两种分析方法计算的相关系数较高而对应的其他方法相关系数较低, 这种波段属于该成分的特征波段。 SOM在450, 580和820 nm附近的波段属于主要特征波段范围; 超过1 400 nm, 例如1 495, 2 210以及2 420 nm左右, 属于SOM的特征波段。 而SMC的特征波段主要分布在大于1 400 nm, 从PCC和MBC曲线可以看出波段均表现出较高的相关性。 1 460和1 950 nm附近是明显的主要特征波段, 与前人总结的SMC特征波段一致。 Fe在整个可见光近红外光谱范围区间, 并没有很明显的主要特征峰, 但较原始反射率相比, 三种分析方法计算得出的所有波段相关度有所提高, PCC的最大值从0.18提升至0.47。 从GRD曲线看, 特征间的相关系数大小差异不明显; 从PCC与MBC曲线判断, 730, 1 360, 1 750以及2 110 nm附近属于Fe的特征波段。 pH的特征波段分散在整个范围内, 从PCC与GRD曲线来看, 特征波段集中在460, 920, 2 240和2 430 nm附近; 从MBC来看, 还集中在1 100和1 740 nm附近。

图5

Fig.5

	Figure Option ViewDownloadNew Window
	图5 SNV-CWT处理前后的反射率与不同成分的关系Fig.5 Correlations between different soil properties and reflectance before and after SNV-CWT

总体分析, 对反射率具有主导影响的属性, 如SOM与SMC, 不同分析方法均可筛选出其对应的主要特征波段。 而影响较弱或含量较低的属性, 如pH与Fe, 主要特征波段不明显, 分析方法不同时所对应的特征波段会有所差别。 从Fe和pH的相关关系图中可以看出, 对于主要特征波段不明显的土壤成分, 不同方法计算的波段相关系数排序不同, 由此所确定的特征波段不同。 因此使用三种分析方法计算波段与土壤成分的相关系数为指标, 采用过滤式方法筛选波段作为模型输入变量并进行模型精度分析。

表3和表4分别记录了以PCC和GRD为过滤式筛选标准的具体建模指标。 从表中可以看出, PCC和GRD对四种成分的计算精度相似; 除SOM外, 其他三种成分所需的波段数目均较多, 且与多维特征建模(表2)相比精度近似。 表中N(Numbers)表示样本数量, T(Threshold)表示相关系数阈值。

表3

Table 3

表3(Table 3)

表3 线性相关法筛选特征及建模精度 Table 3 Statistical indicators of models taking bands selected by PCC as input variables Contents N T Max Min R2 Mse RPD SOM 3 0.80 0.81 0 0.91 0.12 3.49 SMC 99 0.75 0.90 0 0.91 1.66 4.22 Fe 75 0.11 0.22 0 0.43 0.02 1.34 pH 71 0.29 0.44 0 0.65 0.17 1.73

表3 线性相关法筛选特征及建模精度 Table 3 Statistical indicators of models taking bands selected by PCC as input variables

表4

Table 4

表4(Table 4)

表4 灰色关联分析法筛选特征及建模精度 Table 4 Statistical indicators of models taking bands selected by GRA as input variables Contents N T Max Min R2 Mse RPD SOM 2 0.90 0.90 0.68 0.91 0.13 3.32 SMC 89 0.92 0.95 0.74 0.90 1.76 4.02 Fe 80 0.81 0.83 0.74 0.45 0.02 1.38 pH 98 0.82 0.85 0.70 0.65 0.17 1.75

表4 灰色关联分析法筛选特征及建模精度 Table 4 Statistical indicators of models taking bands selected by GRA as input variables

图6为根据MBC值逐步添加特征得到的RF验证模型R²曲线。 与基于PCC和GRD指标分析的特征筛选方法相比, 此方法计算得出的相关关系更适于提高预测模型的精度。 从子图中曲线的变化趋势可以看出, 对于SOM和SMC, 波段数量的增加对模型精度的影响并不大, 基于第6尺度分解系数的SOM与第7尺度分解系数的SMC验证模型R²基本位于0.9~0.95之间, 且可以在较少波段的情况下获得较优的预测模型。

图6

Fig.6

	Figure Option ViewDownloadNew Window
	图6 基于模型的相关关系筛选特征及建模精度Fig.6 Statistical indicators of models taking bands selected by MBC as input variables

Fe的模型精度随着波段数目的增加明显降低。 在可见光近红外波长范围内, 由于SOM含量较高, 在1 400 nm前对反射率具有主导影响, 且SMC对光谱的吸收作用, 干扰反射率曲线, 使1 400 nm后的反射率曲线扭曲, 致使Fe的特征波段信息被掩盖。 同时由于Fe包含的成分复杂, 易受到外界环境及有机质的影响^[5], 较难估算。 相关研究经过100目研磨及MSC处理后的验证模型精度(R²=0.6, RPD=0.86)相比, 图6中Fe的建模当筛选波段数目为8时, 最佳模型R²为0.67, 此时RPD为1.76(1.4< 1.76< 2), 模型较好可用于粗略估算。 pH的预测模型在入选波段数目为41时达到最佳, 此时验证模型R²=0.80, RPD=2.24。