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

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

(吉林大學地球探測科學與技術學院，長春 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，模型可信度較高，驗證了利用隨機森林回歸估算土壤碳酸鈣含量的可行性。