不同粒徑對土壤有機(jī)質(zhì)含量可見-近紅外光譜預(yù)測的影響

鐘翔君, 楊 麗*, 張東興, 崔 濤, 和賢桃, 杜兆輝

1. 中國農(nóng)業(yè)大學(xué)工學(xué)院， 北京 100083

2. 農(nóng)業(yè)農(nóng)村部土壤-機(jī)器-植物系統(tǒng)技術(shù)重點(diǎn)實(shí)驗(yàn)室， 北京 100083

引 言

土壤有機(jī)質(zhì)(soil organic matter， SOM)是表征土壤肥力的關(guān)鍵參數(shù)， 根據(jù)SOM含量進(jìn)行變量播種[1-4]或其他土壤改良作業(yè)對于挖掘土壤潛力、 提高作物產(chǎn)量[5-6]具有重要意義。 傳統(tǒng)以野外采樣及實(shí)驗(yàn)室化學(xué)分析為主的SOM檢測方法精度高， 但工作量大， 成本高， 檢測效率低， 不適于精準(zhǔn)農(nóng)業(yè)大區(qū)域、 快速高效檢測的需求[7]。 近年來， 可見光-近紅外光譜分析技術(shù)因其無損、 簡便、 高效等特點(diǎn)， 在估測SOM含量的研究探索上得到快速發(fā)展， 是土壤特性信息獲取和研究的熱點(diǎn)。

國內(nèi)外許多學(xué)者對SOM與反射光譜之間的關(guān)系[8]、 光譜預(yù)處理方法[9-10]及模型構(gòu)建方法[11-12]開展較多研究， 以實(shí)現(xiàn)SOM含量的高精度快速預(yù)測。 Dhawale等[13]通過對不同質(zhì)地和有機(jī)質(zhì)含量的土樣近紅外光譜進(jìn)行數(shù)學(xué)變換和建模分析， 基于偏最小二乘回歸(partial least squares regression， PLSR)建立了預(yù)測質(zhì)地和有機(jī)質(zhì)含量的光譜校準(zhǔn)方法。 Mouazen等[14]通過對反射率數(shù)據(jù)進(jìn)行不同預(yù)測模型的精度對比， 并確定了該土壤類型下預(yù)測效果最佳的模型。 Dotto等[15]比較了九種模型在七種不同預(yù)處理方法下的SOM預(yù)測精度。 于雷等[16]分析了不同形式的變換光譜與SOM含量的相關(guān)性， 確定連續(xù)統(tǒng)去除和PLSR建立的預(yù)測模型具有顯著的預(yù)測效果。

為了適應(yīng)精準(zhǔn)農(nóng)業(yè)土壤信息高效獲取的要求， 需要針對田間復(fù)雜的作業(yè)環(huán)境開展SOM含量動(dòng)態(tài)檢測的研究， 并開發(fā)適于動(dòng)態(tài)作業(yè)的SOM檢測儀器。 由于土壤的光譜特性是多種土壤屬性交互的綜合反映， 土壤粒徑會(huì)對光譜反射率產(chǎn)生很大影響[17]， 進(jìn)而影響田間SOM含量的動(dòng)態(tài)獲取， 為儀器的開發(fā)增加了難度。 但目前尚未有統(tǒng)一土壤粒徑標(biāo)準(zhǔn)進(jìn)行光譜預(yù)測， 通常為了獲取較好的模型預(yù)測效果， 在測量光譜之前， 先將采集的土壤經(jīng)過不同的處理[18-19]， 如經(jīng)2 mm或1 mm網(wǎng)篩篩分后制備待測樣本， 以降低土壤粒徑對不同土壤屬性預(yù)測的影響， 但不同土壤處理標(biāo)準(zhǔn)得到的預(yù)測結(jié)果間存在差異。 有些研究人員選取過一次篩網(wǎng)后的混合粒徑范圍的土樣(<2 mm或<1 mm)進(jìn)行光譜測量， 不能區(qū)分不同均勻粒徑范圍(如1～2和0.5～1 mm等)對預(yù)測結(jié)果的影響。

針對以上問題， 本工作通過對河北廊坊采集的砂壤土進(jìn)行粒徑分級， 配備1～2， 0.5～1， 0.25～0.5， 0.1～0.25和<0.1 mm五種均勻粒徑及<1 mm的混合粒徑范圍的土樣， 分析SOM含量和光譜反射率的關(guān)系， 并研究不同粒徑土壤光譜預(yù)處理及建模方法對SOM含量的影響， 探究不同粒徑下SOM含量特征波長篩選及提取方法， 并構(gòu)建相應(yīng)的SOM含量預(yù)測模型， 為實(shí)現(xiàn)SOM田間動(dòng)態(tài)預(yù)測及儀器設(shè)計(jì)提供理論及模型參考。

1 實(shí)驗(yàn)部分

1.1 樣本采集

以五點(diǎn)采樣法在河北省廊坊市(39°19′N， 116°17′E)的農(nóng)田區(qū)域中采集0～20 cm的耕層土壤， 該地區(qū)土壤類型主要為砂壤潮土。 采樣時(shí)將表層的石塊及動(dòng)植物殘?bào)w清除， 以保證樣本的完整性。 將采集的樣本充分混合并裝入密封袋中， 作為該區(qū)域的待測土壤樣品， 帶回實(shí)驗(yàn)室進(jìn)行下一步處理。

將田間采集的土樣置于DHG-9123A恒溫干燥箱中， 在40 ℃下烘干至恒重， 將土樣研磨并依次過2， 1， 0.5， 0.25和0.1 mm篩網(wǎng)， 將不同篩網(wǎng)中剩余的土壤封存?zhèn)溆茫?即得到粒徑分別為1～2， 0.5～1， 0.25～0.5， 0.1～0.25和<0.1 mm五種均勻粒徑及<1 mm混合粒徑共計(jì)6種粒徑土樣。 由于<0.1 mm的土樣難以再次進(jìn)行分級， 因此將粒徑<0.1 mm的土樣作為均勻樣品并用于后續(xù)處理。

圖1 不同粒徑的土壤樣本

為獲取不同SOM含量梯度的樣本， 將篩后的土樣分別等量取30份， 并將前期實(shí)驗(yàn)室制備的高有機(jī)質(zhì)含量的營養(yǎng)土(烘干、 研磨、 過篩)分別等比例加入對應(yīng)制備好的粒徑土樣中， 共制備30份不同有機(jī)質(zhì)梯度的土樣， 為儀器設(shè)計(jì)提供完整梯度的SOM含量樣本。 分別將土樣以四分法分為兩份， 主要用于實(shí)驗(yàn)室分析及光譜測試兩部分。

1.2 土壤有機(jī)質(zhì)含量測定

土壤有機(jī)質(zhì)含量采用基于非色散紅外吸收原理的TOC元素分析儀(vario TOC cube， 德國)進(jìn)行測定， 測量時(shí)， 土壤中的有機(jī)物被高溫燃燒氧化后經(jīng)過光學(xué)傳感器， 通過光強(qiáng)的變化測出氣體的濃度， 進(jìn)而實(shí)現(xiàn)總有機(jī)碳的分析， 如圖2(a)所示。 為保證數(shù)據(jù)的有效性， 每個(gè)樣本作5個(gè)重復(fù)， 均值乘以系數(shù)1.724作為最終土樣有機(jī)質(zhì)的有效值， 得到土壤有機(jī)質(zhì)含量統(tǒng)計(jì)結(jié)果如表1所示。

表1 土壤有機(jī)質(zhì)含量統(tǒng)計(jì)

1.3 土壤光譜數(shù)據(jù)采集

采用海洋光學(xué)公司(Ocean Optics, Inc., 美國)的QE Pro高性能光譜儀及NIR Quest512-2.5近紅外光譜儀采集土壤樣品的光譜數(shù)據(jù)， 如圖2(b)所示。 其中， NIR Quest512-2.5近紅外光譜儀光學(xué)分辨率為9.0 nm， 積分時(shí)間為10 ms， 重復(fù)掃描次數(shù)為100次。 QE Pro高性能光纖光譜儀光學(xué)分辨率為1.7 nm， 積分時(shí)間為20 ms， 重復(fù)掃描次數(shù)為50次。 采用全光譜一分三型光纖， 反射光纖分別與對應(yīng)的兩個(gè)光譜儀連接進(jìn)行平行測量， 入射光纖(光纖探頭)與樣品表面的采樣距離為3 mm。 通過筆記本電腦的Ocean View軟件獲取樣本300～2 500 nm波長范圍反射光譜。

分別以四分法取樣的方式， 將不同有機(jī)質(zhì)含量的土壤樣本置于直徑40 mm、 高15 mm的培養(yǎng)皿中， 并用標(biāo)準(zhǔn)塊將表面刮平待測(容重1.32 g·cm-3左右)， 每個(gè)土壤樣本作3個(gè)重復(fù)。 試驗(yàn)前， 用99%漫反射標(biāo)準(zhǔn)白板(Ocean Optics, Inc., 美國)對裝置進(jìn)行校正， 且為降低樣本間環(huán)境誤差， 每采集十個(gè)樣本， 用標(biāo)準(zhǔn)白板校正一次。 為保證樣品光譜數(shù)據(jù)的可靠性， 減小因采樣位置不同造成光譜數(shù)據(jù)的差異， 試驗(yàn)時(shí)， 每個(gè)樣本選取5個(gè)不同位置采集光譜， 每個(gè)位置連續(xù)采集3條光譜， 平均值作為該樣本的反射光譜， 每個(gè)粒徑范圍共得到450條光譜數(shù)據(jù)。

圖2 試驗(yàn)裝置圖

1.4 土壤光譜數(shù)據(jù)預(yù)處理

考慮到光譜儀在低于380 nm和高于2 400 nm波長范圍的光譜信噪比較低， 數(shù)據(jù)噪聲大， 所以剔除這部分?jǐn)?shù)據(jù)， 只保留380～2 400 nm波長范圍作為樣本的有效數(shù)據(jù)。 不同樣本測量時(shí)的實(shí)際操作方式差異不可避免的會(huì)對測量結(jié)果產(chǎn)生影響， 需要對樣本數(shù)據(jù)進(jìn)行去異常值篩選。 通過蒙特卡洛交叉驗(yàn)證法(monte Carlo cross validation， MCCV)對異常樣本進(jìn)行篩選剔除， 得到1～2， 0.5～1， 0.25～0.5， 0.1～0.25， <0.1及<1 mm粒徑土樣的有效光譜數(shù)據(jù)分別為430， 422， 425， 424， 434及426組。

由于環(huán)境和儀器等因素易造成基線漂移， 對所有光譜數(shù)據(jù)采用Savitzky-Golay卷積平滑法進(jìn)行平滑濾波處理， 增強(qiáng)光譜特征， 并將平滑后的光譜作為土樣的原始光譜(R)作后續(xù)處理。 對光譜進(jìn)行相應(yīng)的數(shù)學(xué)變換可將光譜中的有效信息表現(xiàn)出來， 放大相應(yīng)光譜特征。 因此， 對原始光譜(R)分別進(jìn)行倒數(shù)(IR)、 對數(shù)(LR)、 以及一階導(dǎo)數(shù)(FDR)等處理， 并分析不同變換光譜與SOM含量的響應(yīng)關(guān)系。

1.5 基于CARS的SOM特征提取

競爭性自適應(yīng)重加權(quán)算法(competitive adaptive reweighted sampling， CARS)[20]是指在采樣時(shí)利用自適應(yīng)重加權(quán)采樣結(jié)合指數(shù)衰減函數(shù)， 對偏最小二乘回歸(PLSR)模型的回歸系數(shù)進(jìn)行篩選， 去除絕對值較小系數(shù)所在波長， 基于交叉驗(yàn)證計(jì)算并篩選出交叉驗(yàn)證均方根誤差小的波長集合， 該集合包含的變量即為最優(yōu)特征變量組合。 根據(jù)1.4節(jié)中所述， 通過分析不同變換光譜(R， IR， LR， FDR)與SOM含量的響應(yīng)關(guān)系結(jié)果， 選取與SOM含量相關(guān)性較高的光譜進(jìn)行CARS特征變量提取。

1.6 模型建立與檢驗(yàn)

采用光譜-理化值共生距離法(sample set partitioning based on joint x-y distance， SPXY)將70%樣本分為建模集， 30%樣本作為預(yù)測集。 首先， 對不同粒徑的原始光譜(R)及變換光譜(IR， LR， FDR)采用偏最小二乘回歸(PLSR)結(jié)合交叉驗(yàn)證， 建立不同粒徑下不同光譜SOM含量預(yù)測模型， 并分析不同粒徑對SOM含量預(yù)測的影響； 綜合6個(gè)粒徑的原始光譜(R)、 一階導(dǎo)數(shù)(FDR)、 CARS特征變量數(shù)據(jù)作為自變量， SOM含量作為因變量， 分別構(gòu)建基于PLSR的SOM含量預(yù)測模型， 分析不同模型的預(yù)測結(jié)果。 模型的評價(jià)指標(biāo)主要有決定系數(shù)(R2)、 校正均方根誤差(root mean square error of calibration， RMSEC)、 預(yù)測均方根誤差(root mean square error of prediction， RMSEP)及相對分析誤差(residual prediction deviation， RPD)等指標(biāo)。 其中， 決定系數(shù)R2越接近1表明模型的預(yù)測結(jié)果越好， 均方根誤差RMSEP和RMSEC越小越好。 相對分析誤差RPD為預(yù)測集標(biāo)準(zhǔn)差與預(yù)測均方根誤差RMSEP的比值， 當(dāng)RPD≥2.0時(shí)表示模型有較高的可靠性與預(yù)測性能， 1.4

2 結(jié)果與討論

2.1 不同粒徑土壤光譜特征分析

將每種粒徑的土壤光譜數(shù)據(jù)分別按照粒徑及波長進(jìn)行運(yùn)算， 得到土壤粒徑與光譜反射率變化關(guān)系如圖3所示。 從圖3(a)中可以看出， 隨著粒徑的減小， 均勻粒徑平均反射率及變異系數(shù)逐漸增加， 與Sun等[21]的結(jié)論類似。 混合粒徑的平均反射率大于均勻粒徑， 其變異系數(shù)值介于0.25～0.5和0.1～0.25 mm粒徑之間。 從圖3(b)中可以看出， 不同粒徑土樣光譜在不同波長范圍內(nèi)變化趨勢類似， 且隨波長呈現(xiàn)不同變化幅度。 在波長大于540 nm范圍內(nèi)， 平均反射率增幅明顯， 且不同粒徑樣品間的差異也愈加顯著， 在波長為1 410， 1 910及2 210 nm附近位置均出現(xiàn)不同幅度的水分吸收谷[14-15]。 在波長970 nm左右， 由于儀器拼接使得該位置的反射率產(chǎn)生小幅度變化。

圖3 土壤粒徑與光譜反射率變化關(guān)系

2.2 不同粒徑SOM含量與光譜反射率相關(guān)性分析

為分析不同粒徑下光譜變換反射率與SOM含量的相互關(guān)系， 對不同粒徑SOM含量分別與平滑后的R， IR， LR及FDR等變換光譜進(jìn)行相關(guān)性分析， 得到最大相關(guān)系數(shù)與其所在波長位置的結(jié)果如表2所示。 由表2可以看出， 粒徑范圍<1 mm的混合土樣的最大相關(guān)系數(shù)值高于均勻粒徑土樣， 且不同粒徑SOM含量與不同變換光譜反射率的相關(guān)性水平存在顯著差異， 其中FDR變換后SOM含量與光譜反射率的最大相關(guān)系數(shù)明顯大于其他3種光譜變換。

表2 不同粒徑下SOM含量與光譜反射率最大相關(guān)系數(shù)

不同土壤粒徑SOM含量與光譜變換反射率在不同波長的變化關(guān)系如圖4所示。 由于R， IR， LR變換與SOM含量的相關(guān)性變化趨勢類似， 所以只以R、 FDR的結(jié)果為參考進(jìn)行分析。 從圖中可以看出， R變換不同粒徑土樣的相關(guān)系數(shù)在不同波長下的變化趨勢差異不大， 在近紅外波長范圍內(nèi)的相關(guān)性明顯高于可見光波長區(qū)域。 而隨著粒徑的減小， SOM含量與光譜反射率在全波段范圍的相關(guān)性變化也愈加明顯， 其中最大相關(guān)系數(shù)位于<1 mm混合粒徑， 粒徑<0.1 mm的均勻土樣其SOM含量與不同變換光譜相關(guān)系數(shù)的變化率最高。 FDR變換后SOM含量與光譜反射率的相關(guān)性變化明顯， 表明FDR變換可明顯改變?nèi)ǘ畏秶cSOM含量的相關(guān)性[22]， 且粒徑越小， 相關(guān)性系數(shù)變化率越大， 其中均勻粒徑<0.1 mm和<1 mm的混合粒徑的相關(guān)系數(shù)變化范圍明顯大于其他粒徑， 后續(xù)應(yīng)用此光譜進(jìn)行特征提取。 另外， 由于兩臺(tái)光譜儀在970 nm左右進(jìn)行線性拼接， 所以導(dǎo)致在該波長范圍內(nèi)相關(guān)性系數(shù)值產(chǎn)生明顯變化。

2.3 基于CARS的SOM含量特征提取

將6個(gè)粒徑的FDR數(shù)據(jù)分別按照CARS算法進(jìn)行特征變量篩選， 得到不同粒徑下SOM含量的特征變量的分布情況如圖5(a)所示， 從圖中可以看出， 不同粒徑篩選出的特征變量數(shù)及波長位置存在顯著差異， 表明不同粒徑對SOM含量預(yù)測在不同波長范圍影響不同， 也證明了上述粒徑在不同波長范圍對光譜反射率及SOM相關(guān)性影響不同的結(jié)論。 另外， 不同粒徑下的特征變量在某些波長范圍內(nèi)(如390～400和840～850 nm)接近重合， 表明在該波長范圍內(nèi)反射率對SOM含量較為敏感[23]。 為更加直觀描述不同粒徑下特征變量在不同波長的重合次數(shù)， 以波長為橫坐標(biāo)， 以相同波長位置特征變量重合的次數(shù)為縱坐標(biāo)， 繪制如圖5(b)所示分布圖。 從圖中可以看出， 不同粒徑下特征變量在不同波長的重合次數(shù)差異明顯。

圖4 不同粒徑下SOM含量與光譜反射率相關(guān)性

圖5 基于CARS的不同波長SOM含量特征變量分布

依據(jù)上述分析， 將不同粒徑下的特征變量按照變量集中位置與重合次數(shù)進(jìn)行二次提取， 篩選出特征變量在不同粒徑下分布較為集中的波長作為SOM含量預(yù)測的特征波長。 以圖5(b)得到的在相同波長下不同粒徑特征變量重合次數(shù)為4和5的波長為中心， 配合圖5(a)中不同粒徑特征變量分布較為集中的區(qū)域， 同時(shí)對不同粒徑下特征變量分布較為集中的重合2和3次的波長進(jìn)行篩選， 共篩選出162組特征波長， 特征波長數(shù)量降低至全波段數(shù)量的13.1%， 篩選出的特征波長分布結(jié)果如圖6所示， 從圖中可以看出， 特征波長大多集中于近紅外區(qū)域[24]。

2.4 不同粒徑SOM含量預(yù)測模型構(gòu)建與分析

采用SPXY法對不同粒徑的樣本進(jìn)行劃分， 將70%樣本分為建模集， 30%樣本作為預(yù)測集。 首先， 分別以不同粒徑梯度下R， IR， LR及FDR等4種變換光譜反射率為自變量， SOM含量為因變量， 采用PLSR建立不同粒徑下的SOM預(yù)測模型， 分析不同粒徑及不同光譜變換形式對SOM含量預(yù)測的影響； 其次， 以篩選出的CARS特征變量數(shù)據(jù)以及全部粒徑的原始光譜(R)、 一階導(dǎo)數(shù)(FDR)作為自變量， SOM含量作為因變量， 構(gòu)建PLSR預(yù)測模型。 模型結(jié)果如圖7和圖8所示。

圖6 基于CARS篩選的SOM特征分布

圖7 不同粒徑下SOM含量PLSR預(yù)測模型構(gòu)建與評價(jià)

由于土壤光譜反射率受多種因素(如質(zhì)地、 顏色等)的影響， 本研究中未針對這些影響因素開展研究， 在后續(xù)研究中， 需要考慮影響土壤光譜的其他因素及不同因素間的交互對建模精度的影響。 另外， 為了適應(yīng)精準(zhǔn)農(nóng)業(yè)土壤信息高效獲取的要求， 需要針對田間復(fù)雜的作業(yè)環(huán)境(光照、 作業(yè)速度等)開展SOM含量田間動(dòng)態(tài)檢測的研究， 對數(shù)據(jù)處理及建模方法進(jìn)一步優(yōu)化， 以實(shí)現(xiàn)田間SOM動(dòng)態(tài)高精度檢測。

3 結(jié) 論

通過對田間采集的砂壤土進(jìn)行粒徑分析， 分別制備了1～2， 0.5～1， 0.25～0.5， 0.1～0.25和<0.1 mm五種均勻粒徑及<1 mm的混合粒徑范圍的土壤， 比較了不同粒徑樣品的光譜反射率差異， 并對原始光譜R進(jìn)行倒數(shù)IR、 對數(shù)LR、 一階導(dǎo)數(shù)FDR等3種光譜變換， 分析了不同粒徑SOM含量與不同變換光譜反射率的相關(guān)性， 基于CARS方法選取SOM含量的特征變量， 結(jié)合PLSR分別建立了相應(yīng)的SOM含量預(yù)測模型。 主要結(jié)論如下：

(1)不同粒徑土樣平均光譜反射率隨著粒徑的減小逐漸增加， 其中粒徑<0.1 mm的均勻粒徑和<1 mm混合粒徑的土樣的平均反射率最大； 在不同波長范圍內(nèi)不同粒徑土樣的平均反射率呈現(xiàn)不同的變化幅度， 尤其波長大于540 nm范圍， 粒徑間的差異更加明顯。

(2)隨著粒徑的減小， SOM含量與光譜反射率在全波段范圍的相關(guān)性變化也愈加明顯， 且FDR變換可明顯改變?nèi)ǘ畏秶cSOM含量的相關(guān)性， 其中粒徑<1 mm和<0.1 mm土樣的相關(guān)系數(shù)變化幅度最為明顯。

圖8 不同PLSR模型預(yù)測值與實(shí)測值的比較