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

吳倩，姜琦剛，史鵬飛，張莉莉

(吉林大學(xué)地球探測科學(xué)與技術(shù)學(xué)院，長春 130026)

0 引言

土壤是自然界不可或缺的資源[1]，它為作物的生長提供了一定的養(yǎng)分和水分。碳酸鈣(CaCO3)是土壤中碳酸鹽的主要成分,對土壤物理、化學(xué)性質(zhì)、碳的沉積、土壤酸堿度、土壤膠體性狀等的研究都有重要意義[2-3]。近年來，隨著遙感科學(xué)的發(fā)展，高光譜遙感應(yīng)用的深入，人們越來越體會到定量遙感的必要性，基于高光譜的土壤理化性質(zhì)快速、定量、準(zhǔn)確的模型構(gòu)建成為農(nóng)業(yè)遙感研究領(lǐng)域研究的重點[4]。土壤的光譜反射特性是土壤基本特性之一，與土壤中的各種理化性質(zhì)有著十分密切的關(guān)系，土壤中各種物質(zhì)特殊的光譜特征以及土壤組成物質(zhì)的多樣性都會在光譜曲線上體現(xiàn)，碳酸鈣的含量對土壤光譜反射特性具有一定的影響，Gaffey[5]研究發(fā)現(xiàn)碳酸鹽礦物由于碳酸根的振動過程，在可見-近紅外(350～2 550 nm )光譜區(qū)域有7個較強的吸收特征； 張婷婷[6]研究認(rèn)為土壤中碳酸鹽含量對土壤的光譜特征有一定的影響，在2 230～2 270 nm之間有明顯的吸收特征； Stenberg[7]研究認(rèn)為土壤碳酸鹽在2 335 nm處有較強的吸收特征； 以往研究表明碳酸鈣的可見光、近紅外(VNIR，400～1 200 nm)、短波紅外(SWIR，1 200～2 500 nm)光譜可作為土壤碳酸鈣含量估算的理論依據(jù)，是一種物理性的無損傷和低成本的研究方法，可以用來區(qū)分土壤類型和評估大量土壤特性。為了更好地理解電磁波與地表特征之間的相互作用，國內(nèi)外學(xué)者將經(jīng)典的數(shù)學(xué)物理理論與遙感實踐相結(jié)合，建立了近百種定量遙感模型，這些模型大體可分為統(tǒng)計模型、物理模型、半經(jīng)驗?zāi)Ｐ?種。利用高光譜數(shù)據(jù)估測土壤理化性質(zhì)的模型有很多，近年來的研究中大多運用統(tǒng)計模型，建立遙感數(shù)據(jù)與地面觀測數(shù)據(jù)之間的線性回歸方程，統(tǒng)計模型簡便易行，一般僅包含3～6個參數(shù)。以往研究中多運用偏最小二乘回歸、多元線性同步回歸、支持向量機等方法，且研究內(nèi)容大多集中在對土壤有機質(zhì)含量、土壤水分、土壤游離鐵含量的估算建模上。本次研究選取隨機森林回歸算法(random forest regression，RFR)對土壤碳酸鈣含量進行估算，隨機森林回歸是一種集成學(xué)習(xí)模型，其優(yōu)化了組合機器學(xué)習(xí)和CART決策樹算法，并且對異常值和有噪聲的數(shù)據(jù)具有一定的容忍度[8]，目前更多的應(yīng)用于數(shù)據(jù)分類和回歸分析[9]。

研究選取的土壤類型為黃綿土，黃綿土質(zhì)地疏松，呈淺棕色，剖面發(fā)育不明顯。廣泛分布于甘肅、陜西、山西等水土流失較嚴(yán)重地區(qū)，是陜西省分布面積最大的土壤類型，其肥力低下，主要分布在坡度大、植被較少地區(qū)，是造成土壤養(yǎng)分有效性低、pH高，影響土壤結(jié)構(gòu)及加速土壤侵蝕的重要因素之一，因此基于高光譜的土壤碳酸鈣含量估算可以準(zhǔn)確快速地獲取土壤屬性信息以及為水土流失地區(qū)提供治理等的理論參考，對于評價成土因素以及對土壤性質(zhì)的研究具有重要意義。

1 材料與方法

1.1 研究區(qū)概況

研究區(qū)主要位于陜西省北部的黃土高原地區(qū)。107°15′41″～111°14′40″E，34°12′21″～39°29′17″N，黃土高原屬于干旱區(qū)與半干旱區(qū)； 在氣候分區(qū)上，該地區(qū)橫跨溫帶和暖溫帶，年平均降水量約在300～700 mm之間，長期的水土流失導(dǎo)致生態(tài)環(huán)境較脆弱。陜西省土壤類型眾多，本次研究的土壤類型為黃土高原上分布面積最大的黃綿土，研究共采集了78個土樣，其中28個位于陜西省最北部的榆林市，41個位于延安市，9個位于渭南市。圖1 為研究區(qū)地理位置及采樣點分布。

圖1 研究區(qū)地理位置及采樣點分布Fig.1 Location of the study area and the distribution of sampling sites

1.2 樣本采集與分析

根據(jù)陜西省黃綿土空間分布特征，于2016年7—8月期間在榆林市、延安市、渭南市共采集具有典型代表性的78個黃綿土土壤樣品，對土壤樣本進行風(fēng)干、碾磨及過1 mm孔篩等操作后，采用四分法，將每個樣品分為2份，一份利用氣量法進行實驗室土壤碳酸鈣含量測定，另一份用作土壤反射光譜測定。

1.3 光譜測定

土壤VNIR/SWIR反射光譜是在暗室中控制輻照度條件測量的，測定采用了美國 Spectra Vista 公司生產(chǎn)的SVC HR-1024i便攜式光譜儀進行，其波長范圍為350～2 500 nm。在0.6～1 m之間放置1盞50 W白熾燈，設(shè)置60°方位角和15°天頂角，這種配置確保了由于陰影造成的最小擾動，使測量結(jié)果不那么依賴于粗糙度。用8°透鏡提供直徑為0.13 m的光斑。測量中利用白板進行校正，獲得了3個連續(xù)的和50個單獨反射測量系列，噪聲水平低于0.001，光源的照射方向與垂直照射之間的夾角為30°，測量土樣在4個方向上的光譜曲線。

1.4 數(shù)據(jù)預(yù)處理

為了提高光譜數(shù)據(jù)的質(zhì)量，研究利用ENVI5.3對建模前的數(shù)據(jù)進行光譜數(shù)據(jù)庫的建立以及10 nm間隔的重采樣，重采樣處理后每個土壤樣本的光譜數(shù)據(jù)曲線包含216個波段。重采樣后的光譜曲線依然存在一定量噪聲，為了提高預(yù)測模型的精度，采用Savitzky-Golay 卷積平滑法(S-G)對數(shù)據(jù)進行進一步的處理[10]，平滑去噪后，對土壤碳酸鈣光譜數(shù)據(jù)進行數(shù)學(xué)變換，得到原始光譜(R)的一階微分(First order differential，R″)數(shù)據(jù)，二階微分(Second order differential，R″)數(shù)據(jù)，連續(xù)統(tǒng)去除(Continuum Removal，CR)數(shù)據(jù)。

1.5 敏感波段選取

采用相關(guān)分析法(analysis of correlation，CA)與連續(xù)投影算法(successive projections algorithm，SPA)分別進行土壤碳酸鈣含量與光譜反射率的相關(guān)性分析以及敏感波段的篩選，并依據(jù)建模結(jié)果對兩種方法進行評價。

CA是指對兩個或多個具備相關(guān)性的變量元素進行分析，從而衡量兩個變量因素的相關(guān)密切程度[10]。 需要將4種光譜處理方法與所測樣本的碳酸鈣含量的真實值逐個波長進行相關(guān)性分析，并進行p=0.01顯著性檢驗，根據(jù)相關(guān)性的高低來確定敏感波段，相關(guān)性越高，波段越敏感。相關(guān)分析采用SPSS中的相關(guān)性函數(shù)進行，通過顯著性檢驗的波段作為接下來的建模的敏感波段。相關(guān)系數(shù)公式如下：

(1)

SPA在多自變量校正領(lǐng)域是一種采用前向選取的特征波段變量選擇的方法[12]。能夠有效地減少信息重疊，使變量之間的共線性達到最小，大大減少了建模變量的個數(shù)，有效提高建模的效率。實現(xiàn)用較少的信息量代表多數(shù)樣本的光譜信息的效果[10]，SPA以均方根誤差(ROOT mean square error，RMSE)為評價指標(biāo)，將RMSE最小值下的波長個數(shù)確定為敏感波段數(shù)[13-14]，連續(xù)投影算法在Matlab2018a中實現(xiàn)。

1.6 模型建立與精度驗證方法

隨機森林(RFR)是一種集成算法(Ensemble Learning)，這種算法優(yōu)化了組合機器學(xué)習(xí)和CART決策樹算法[15]。隨機森林有很多的優(yōu)點，首先它能在決策類別時，去評估變量的重要性，其次它對于數(shù)據(jù)集較多的輸入變量的處理精度較好，隨機森林的兩個應(yīng)用主要是分類與回歸。近幾年隨機森林在回歸分析方面的應(yīng)用也越來越廣泛。該算法的估算精度比較高，樹形結(jié)構(gòu)自由生長，能夠有效避免其他估算模型中的過度擬合現(xiàn)象，是對傳統(tǒng)的決策樹模型的一種優(yōu)化，對數(shù)據(jù)量較大的數(shù)據(jù)集比較適用[15]。隨機森林回歸模型的構(gòu)建運用R.3.4.1中的Random Forest包進行。

研究選取土壤碳酸鈣含量實測值與估測值的均方根誤差(RMSE)，決定系數(shù)(R2)以及相對分析誤差(residual predictive deviation，RPD)3個參數(shù)對模型進行精度驗證，均方根誤差與相對分析誤差越小，決定系數(shù)越大，模型的精度越高。

2 結(jié)果與分析

2.1 土壤碳酸鈣描述性統(tǒng)計

首先對全部78個土壤碳酸鈣樣本進行異常值分析，將樣本中碳酸鈣含量不在±3倍標(biāo)準(zhǔn)差中的極度異常的4個樣本剔除，從剩余的74個樣本中選擇出52個樣本用于建模，22個樣本用于驗證，表1為土壤樣本碳酸鈣含量的統(tǒng)計結(jié)果，樣本總體碳酸鈣含量范圍為11.109～175.077 g·kg-1，其中建模集中的52個碳酸鈣含量范圍為11.109～173.509 g·kg-1，驗證樣本碳酸鈣含量范圍為11.236～175.077 g·kg-1，樣本的選取基本上涵蓋了總體樣本中的數(shù)據(jù)范圍。建模集與驗證集的偏度為0.46和0.25，均在0值附近，基本服從正態(tài)分布，建模集、驗證集與總體樣本的偏度差異不大，樣本適合用于建立模型。

表1 土壤樣本碳酸鈣含量的統(tǒng)計結(jié)果Tab.1 Statistics of calcium carbonate in soil samples

2.2 土壤碳酸鈣含量與光譜相關(guān)性分析

圖2(a)為全部碳酸鈣樣本波長與反射率之間的原始光譜圖，由圖中可以看出，不同碳酸鈣含量的土壤反射率曲線總體趨勢基本相同，深度有所不同，吸收特征基本相同，在1 400 nm和1 900 nm和2 200 nm處有明顯的吸收特征，這與前人分析結(jié)果相似，通常認(rèn)為是由OH，AL-OH等引起的[16-17]，除此之外，近紅外波段的反射率總體大于可見光部分的反射率； 如圖2(b)所示，隨著碳酸鈣含量的增加，光譜反射率也隨之增加，兩者呈現(xiàn)正相關(guān)態(tài)勢，且近紅外波段的增長速度明顯高于可見光部分的增長速度，反射率曲線在于2 350 nm處有一吸收谷，這主要是由于碳酸根的影響作用引起的[16]。

圖3為基于相關(guān)分析算法的敏感選波段篩選圖，由圖3可以看出，原始光譜曲線在全波段的相關(guān)性較差，相關(guān)系數(shù)最高值在2 340 nm左右，相關(guān)系數(shù)僅為0.37，曲線整體呈現(xiàn)正相關(guān)的態(tài)勢。經(jīng)過一階微分、二階微分、連續(xù)統(tǒng)去除后，光譜曲線與土壤碳酸鈣含量的相關(guān)性得到了顯著的提高，從而可以進一步提高模型的預(yù)測能力。其中一階微分相關(guān)系數(shù)最大值在1 610 nm處達到了-0.66，二階微分相關(guān)系數(shù)在2 230 nm處達到了0.56，微分運算雖然可以較好地消除背景信號，但高階微分在應(yīng)用中往往表現(xiàn)出不穩(wěn)定性，在應(yīng)用中存在一定的問題，以往研究發(fā)現(xiàn)，在光譜分析中隨著階數(shù)的增加，相關(guān)系數(shù)呈現(xiàn)先增加后減少的趨勢，并在二階處達到最大值[18]，因此，研究采用一階與二階微分。連續(xù)統(tǒng)去除相關(guān)系數(shù)在2 150 nm處為到0.53，均達到中度相關(guān)水平，4種數(shù)學(xué)處理方法分別有54，82，54，44個波段入選建模波段。

圖4為基于連續(xù)投影算法的敏感波段篩選圖，經(jīng)過連續(xù)投影算法降維處理過后篩選的敏感波段在圖中用方塊圈出，經(jīng)S-G平滑后的原始光譜數(shù)據(jù)，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個波段入選； 經(jīng)一階微分后的光譜數(shù)據(jù)，RMSE最小值為26.194 5，共420 nm，470 nm，620 nm，1 510 nm和1 530 nm等9個波段入選； 經(jīng)二階微分后的光譜數(shù)據(jù)，RMSE最小值為21.272 4，共360 nm，420 nm，830 nm，880 nm和1 020 nm等16個波段入選； 經(jīng)連續(xù)統(tǒng)去除后的光譜數(shù)據(jù)，RMSE最小值為22.853 3，共410～430 nm，540 nm和570 nm等15個波段入選。與相關(guān)分析法相比大大減少了模型的輸入變量，可有效提高模型的計算速度。

2.3 建模結(jié)果分析

表2為土壤碳酸鈣含量敏感波段RFR模型，研究表明，經(jīng)數(shù)學(xué)變換的光譜為變量建立的模型精度明顯高于原始光譜，數(shù)學(xué)變換有效提高了土壤碳酸鈣含量與光譜反射率之間的相關(guān)性，驗證集R2均大于0.58，其中二階微分的估算精度最高，驗證集R2為0.82，RPD為2.37，表明模型較為適用； 從敏感波段選取方法看，除基于CR的建模集精度相關(guān)分析法較高于連續(xù)投影算法外，其余數(shù)學(xué)變換方法無論建模集驗證集基于連續(xù)投影算法建模精度均高于相關(guān)分析法，由此看出，連續(xù)投影算法不僅變量少效率高，并且在波段篩選方面較為適用； 就模型本身而言，隨機森林回歸整體的預(yù)測精度較好，驗證了利用隨機森林回歸估算土壤碳酸鈣含量的可行性。就模型的穩(wěn)定性而言，建模集整體精度均高于驗證集，模型缺乏一定的穩(wěn)定性，可能原因為： ①土壤碳酸鈣樣品普遍存在標(biāo)準(zhǔn)差較大、變異性較大的現(xiàn)象，碳酸鈣的含量較為分散，在對碳酸鈣含量進行估算研究中數(shù)據(jù)普遍存在以上問題； ②獲取敏感波段的方式不同。例如，洪長喬等[17]除了利用相關(guān)分析法選取的波段，還添加了前人研究結(jié)果較好的幾個數(shù)據(jù)集作為自己研究的建模集之一，并分別建模進行比較； Gomez等[19]在研究拉比因河谷土壤碳酸鈣含量時，利用可變重要性(VIP分?jǐn)?shù))和B系數(shù)的值來確定光譜帶，由此可看出敏感波段的選取方式是多樣的，選擇方式的不同對模型精度與穩(wěn)定性會產(chǎn)生一定的影響。圖5為基于SPA的4種數(shù)學(xué)變換RFR模型實測值與預(yù)測值擬合圖，可以明顯看出二階微分的模型擬合效果較好。

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

3 結(jié)論

研究采用隨機森林回歸的方法建立了陜西省北部黃綿土高光譜碳酸鈣估算模型，旨在探討不同光譜處理方法與不同敏感波段篩選方法模型構(gòu)建的優(yōu)劣以及利用隨機森林回歸估算碳酸鈣含量的可行性，主要取得了以下3點結(jié)論：

1)碳酸鈣較明顯的吸收特征位于1 610 nm，2 220 nm，2 340 nm，2 270 nm等處，且土壤碳酸鈣含量與全波段光譜呈現(xiàn)較明顯的正相關(guān)關(guān)系。

2)基于連續(xù)投影算法的模型構(gòu)建精度除建模集CR變換外精度均高于相關(guān)分析法，表明連續(xù)投影算法在波段篩選中較為適用。

3)基于連續(xù)投影算法與二階微分變換的隨機森林回歸建模精度最高，建模集R2為0.89，驗證集R2值為0.82，且RPD值大于2，模型可信度較高，驗證了利用隨機森林回歸估算土壤碳酸鈣含量的可行性。