土壤有機質(zhì)含量地面高光譜估測模型對比分析

王永敏， 李西燦， 田林亞， 賈斌， 楊惠

(1.河海大學(xué)地球科學(xué)與工程學(xué)院，南京 211100； 2.山東農(nóng)業(yè)大學(xué)信息科學(xué)與工程學(xué)院，泰安 271018； 3.中建四局第三建筑工程有限公司，遵義 563000； 4.東南大學(xué)交通學(xué)院，南京 210000)

0 引言

土壤有機質(zhì)含量是農(nóng)作物生長的必要物質(zhì)之一，對其進行快速、準確測定和估計對發(fā)展精細農(nóng)業(yè)和提高農(nóng)田管理的科學(xué)性具有重要意義[1]。傳統(tǒng)的土壤有機質(zhì)含量測定方式手段繁瑣，耗時耗力，且有一定的局限性。近年來，高光譜技術(shù)的發(fā)展，使得快速實現(xiàn)土壤有機質(zhì)含量估測成為可能，因而研究合理有效的土壤有機質(zhì)含量高光譜估測方法非常重要。邱華旭等[2]通過對光譜反射率做常規(guī)模型變換，發(fā)現(xiàn)適當(dāng)?shù)哪Ｐ妥儞Q可以使土壤性質(zhì)指標(biāo)與反射率數(shù)據(jù)的相關(guān)性得到有效提高； 郭超凡等[3]和Kokaly等[4]對光譜吸收特性采用連續(xù)移除法和多元逐步線性回歸法處理，可以提高高光譜估測模型的精度； 陳增文等[5]發(fā)現(xiàn)土壤有機質(zhì)含量與光譜反射率呈負相關(guān)，為利用土壤高光譜反射率數(shù)據(jù)進行有機質(zhì)含量估測提供了理論支持和方法引導(dǎo)。但是，諸多研究在光譜敏感波段方面的結(jié)論并不一致，彭杰等[6]研究土壤理化特征發(fā)現(xiàn)，350～2 500 nm波段可以反映土壤特性的細微差異，可用于反演土壤理化參數(shù)； 廖欽洪等[7]和Liu等[8]采用一階微分歸一化技術(shù)處理黑土的反射率數(shù)據(jù)，研究發(fā)現(xiàn)光譜波段范圍為620～810 nm時對有機質(zhì)敏感； Karnieli等[9]發(fā)現(xiàn)在1 720 nm,2 180 nm和2 309 nm處光譜曲線存在吸收特性。此外，采用不同的數(shù)據(jù)處理方法，所建模型精度也存在差異，沈潤平等[10]通過人工神經(jīng)網(wǎng)絡(luò)方法建立土壤有機質(zhì)含量估測模型，總均方根誤差達到了1.31； 于雷等[11]通過連續(xù)小波變換和支持向量機回歸建立的估測模型決定系數(shù)達到了0.83； 彭小婷等[12]采用包絡(luò)線去除法，結(jié)合偏最小二乘回歸法，建立模型的決定系數(shù)達到了0.763。但是，目前綜合比較各種土壤有機質(zhì)含量高光譜估測方法的研究較少。

為探究不同方法用于土壤有機質(zhì)含量高光譜估測的精度，確定有效的估測建模方法，本文通過小波分析和包絡(luò)線去除法實現(xiàn)去噪和參數(shù)提取，采用9種數(shù)據(jù)變換方法，對提取的建模因子進行變換并結(jié)合神經(jīng)網(wǎng)絡(luò)法、多元線性回歸法和最小二乘回歸法共建立了27種地面高光譜土壤有機質(zhì)含量估測模型，比較不同方法所建模型精度，最終確定有效的估測建模方法，以期為相關(guān)研究提供一定的參考。

1 數(shù)據(jù)與方法

1.1 數(shù)據(jù)采集

本文所用數(shù)據(jù)為陜西省橫山縣土壤樣本的地面實測高光譜反射率和有機質(zhì)含量數(shù)據(jù)。陜西省橫山縣(E108°56′～110°02′，N37°22′～38°14′)位于陜西省北部，其地勢西南高，東北低，中間高，南部低，呈傾斜趨勢，海拔約890～1 534 m； 由于處在鄂爾多斯草原與黃土高原的過渡地帶，水土流失嚴重，“三沙七丘”地貌景觀獨特，主要分布有黃綿土、綿沙土、黑壚土和風(fēng)沙土等土壤類型。

土樣采集利用梅花樁取樣法，保證采樣點盡可能均勻分布在研究區(qū)內(nèi)，顧及不同土壤類型，在每個測區(qū)內(nèi)隨機選擇4～5個具有代表性的測點(白界鄉(xiāng)、岔龍灣、高鎮(zhèn)和石灣鎮(zhèn)4個采樣區(qū)測點數(shù)為5，其余為4)，選擇天氣晴朗、光照穩(wěn)定、土壤干燥時除掉測點上方的腐殖落葉層并進行5次光譜測定取平均值作為樣本點的光譜反射率值，然后取15 cm左右深度的表層土，共采集84個樣本，樣本采集區(qū)分布情況如圖1所示。

圖1 土壤樣本采集區(qū)域分布示意圖Fig.1 Area of collecting soil samples

光譜反射率采集使用ASD FieldSpec FR野外光譜儀，采用3°視場角，面對太陽方向，探頭距離目標(biāo)20 cm進行垂直觀測，觀測范圍直徑近2 cm。光譜反射率波段范圍為350～2 500 nm，其中350～1 000 nm光譜分辨率為3 nm，1 000～2 500 nm光譜分辨率為10 nm。對土壤樣本的有機質(zhì)含量和水含量在室內(nèi)進行分析測定。土壤樣本的理化參數(shù)統(tǒng)計結(jié)果如表1所示。在數(shù)據(jù)處理之前，通過比較分析，剔除了數(shù)據(jù)噪聲較嚴重的4個樣本，剩余80個樣本數(shù)據(jù)用于本研究。經(jīng)多次實驗，隨機選擇50個樣本用于建模，剩余30個樣本用于檢驗。

表1 土壤樣本有機質(zhì)含量和水含量統(tǒng)計Tab.1 Statistical characteristics of soil organic matter content and water content (%)

1.2 數(shù)據(jù)處理

1.2.1 小波分析去噪

在數(shù)據(jù)獲取過程中，由于測定方法、外界環(huán)境和操作等原因，獲得的光譜反射率數(shù)據(jù)不可避免地存在噪聲，因此對原始數(shù)據(jù)進行了去噪處理。小波分析是常用且效果較好的去噪方法之一，本文通過Matlab編程實現(xiàn)光譜反射率的小波去噪。選用Bior1.3小波基函數(shù)進行3層小波分解，當(dāng)小波系數(shù)小于閾值時，既不采用硬閾值法置零，也不采用軟閾值法按照固定向量向零附近收縮，而是通過歸一化相關(guān)系數(shù)設(shè)置默認閾值[13]，歸一化相關(guān)系數(shù)公式為

(1)

(2)

(3)

cwJ,K=wJ,K×wJ+1,K，

(4)

式中：cwJ,K為尺度J上K點處的相關(guān)系數(shù)；wJ,K和wJ+1,K分別為小波系數(shù)，由信號小波系數(shù)和噪聲小波系數(shù)2部分組成。通過默認閾值對信號進行重構(gòu)，實現(xiàn)原始光譜反射率數(shù)據(jù)的去噪。

1.2.2 包絡(luò)線去除法提取參數(shù)

由于原始的光譜反射率數(shù)據(jù)量大，本文采用包絡(luò)線去除法壓縮數(shù)據(jù)量并提取參數(shù)。利用VB編程并結(jié)合Excel表格，選擇與有機質(zhì)含量相關(guān)性高的連續(xù)波段為敏感波段，并在選取的敏感波段上進行包絡(luò)線去除，共提取了高度面積、深度面積、高差面積、最大高度、最大深度、最大高差、對稱度和位置等參數(shù)。包絡(luò)線提取參數(shù)的運算公式分別為

C=Ri/[Rstart+T×(λi-λstart)]，

(5)

K=(Rend-Rstart)/(λend-λstart)，

(6)

H=(1/C)-1，

(7)

P=1-C，

(8)

D=(1/C)-C，

(9)

(10)

(11)

(12)

(13)

(14)

(15)

DS=K/(AP-K)，

(16)

(17)

式中：C為包絡(luò)線結(jié)果值；Ri為i處的原始光譜反射率數(shù)據(jù)；Rstart和Rend分別為起止點處的光譜反射率；λi為i處的波長；λstart和λend分別為起止點處的波長；T為待處理的光譜反射率曲線起止點間的斜率；K為最大深度處深度值；H，P和D分別為高度、深度和高差；AH為高度面積；AP為深度面積；AD為高差面積；MH為最大高度；MP為最大深度；MD為最大高差；DS為對稱度；LM為位置；s和e分別為起點和終點位置。按照相關(guān)性由高到低的順序，選取其中相關(guān)性較高的6個特征作為建模因子。

1.2.3 數(shù)據(jù)變換

為提高建模因子與土壤有機質(zhì)含量的相關(guān)性，選擇9種數(shù)據(jù)變換方法進行變換處理，包括平方變換(square transform，SQ)、倒數(shù)變換(reciprocal transform，RE)、指數(shù)變換(exponential transform，EXP)、對數(shù)變換(log transform，LOG)、指數(shù)的平方變換(exponential square transform，ES)、對數(shù)的平方變換(logarithm square transform，LS)、一階微分變換(first order differential transform，DE1)、對數(shù)一階微分變換(logarithm first order differential transform，LGD1)、指數(shù)一階微分變換(exponential first order differential transform，EXD1)，其中一階微分用差分表示。由于位置因子的指數(shù)變換數(shù)值過大，先將位置因子縮小100倍后進行指數(shù)變換。

1.3 建模及精度評定方法

采用BP(back propagation)神經(jīng)網(wǎng)絡(luò)法、多元線性回歸法(multivariable linear regression，MLR)和最小二乘回歸法(least squares regression，LSR)，結(jié)合9種數(shù)據(jù)變換方法，共建立了27種土壤有機質(zhì)含量地面高光譜估測模型。

網(wǎng)絡(luò)參數(shù)設(shè)置對神經(jīng)網(wǎng)絡(luò)模型的結(jié)果質(zhì)量有很大影響，因此選擇合適的參數(shù)(隱含層數(shù)、學(xué)習(xí)速率、平滑因子、學(xué)習(xí)誤差和分級迭代級數(shù))對模型的精度很重要。由于目前沒有很成熟的參數(shù)設(shè)置模型可供使用，本文通過大量實驗，反復(fù)對比，最終確定本實驗的參數(shù)設(shè)置為： 輸入層數(shù)為6，對應(yīng)提取的6個建模因子所對應(yīng)的數(shù)據(jù)； 輸出層數(shù)為1，對應(yīng)樣本有機質(zhì)含量，人為確定前50個樣本數(shù)據(jù)作為學(xué)習(xí)樣本，后30個樣本數(shù)據(jù)作為檢驗樣本； 隱含層數(shù)為15，學(xué)習(xí)速率為1.5 h，平滑因子為0.7，學(xué)習(xí)誤差為0.03，分級迭代級數(shù)為14。通過編程實現(xiàn)MLR和LSR建模。

采用模型的決定系數(shù)R2、檢驗樣本的均方根誤差(root mean square，RMSE)和估測模型的差異顯著性sig值作為估測模型的精度和可靠性評定因子。各模型的sig值由SPSS軟件獲得，用以檢驗?zāi)Ｐ偷娘@著性； 根據(jù)R2和RMSE，比較不同模型的精度和可靠性，最終選出最佳模型。其中，sig<0.05則說明模型是顯著的，R2越大，模型的可靠性越強、模型越穩(wěn)定；RMSE越小，模型精度越高。

2 結(jié)果與分析

2.1 小波分析去噪結(jié)果

通過Bior1.3小波基函數(shù)對原始光譜反射率數(shù)據(jù)進行3層小波分解，實現(xiàn)去噪處理(圖2)。

(a) 小波去噪前光譜反射率(b) 小波去噪后光譜反射率(c) 去噪前后反射率與有機質(zhì)含量相關(guān)系數(shù)

圖2小波去噪前后光譜反射率與有機質(zhì)含量相關(guān)系數(shù)對比

Fig.2Comparisonofcorrelationcoefficientbetweenspectralreflectanceoforiginalandwaveletdenoisingandorganicmattercontent

圖2(a)和(b)中可以看出(以77號樣本為例，其他樣本處理效果與此類似)，小波處理使得原始數(shù)據(jù)中的噪聲得到去除，在1 400 nm附近和2 000～2 500 nm表現(xiàn)明顯。由圖2(c)可以發(fā)現(xiàn)原始光譜反射率數(shù)據(jù)與有機質(zhì)含量呈負相關(guān)且相關(guān)性不高，小波分析去噪提高了光譜反射率與有機質(zhì)含量的相關(guān)性。因此，采用小波去噪后的光譜反射率數(shù)據(jù)代替原始數(shù)據(jù)建立估測模型。

2.2 參數(shù)提取

由圖2的數(shù)據(jù)得知，在波段590～770 nm之間，光譜反射率與有機質(zhì)含量相關(guān)系數(shù)較高，是光譜敏感波段，為減小數(shù)據(jù)處理量，選擇波段600～700 nm之間的光譜反射率數(shù)據(jù)，編程實現(xiàn)包絡(luò)線去除法，并對參數(shù)進行提取。以66號土壤樣本為例，圖3為包絡(luò)線去除前后光譜反射率對比，其他樣本處理后的曲線趨勢與此類似。

圖3 包絡(luò)線去除結(jié)果對比Fig.3 Continuum removal results

實驗共提取高度面積、深度面積、高差面積、最大高度、最大深度、最大高差、對稱度和位置等多項指標(biāo)因子。按照相關(guān)性由高到低的順序，選取位置、對稱度、深度面積、最大深度、最大高差和最大高度等6個參數(shù)作為建模因子，但是各個因子與有機質(zhì)含量相關(guān)性不高，不能直接用于建模，故采用不同數(shù)據(jù)變換方法對各因子進行變換處理(表2)。

表2 相關(guān)系數(shù)對比Tab.2 Comparison of correlation coefficient

2.3 建模因子變換結(jié)果比較

由于原始因子與有機質(zhì)含量的相關(guān)性不高，通過1.2.3節(jié)中的方法對建模因子數(shù)據(jù)進行變換處理。計算各因子與有機質(zhì)含量的相關(guān)系數(shù)(表2)結(jié)果表明，提取參數(shù)經(jīng)LS變換后，除最大深度外，其他參數(shù)與有機質(zhì)含量的相關(guān)系數(shù)均達到0.5以上，能有效提高提取參數(shù)與土壤有機質(zhì)含量的相關(guān)性。此外，EXP，ES，LGD1和EXD1的變換效果也較好。

2.4 建模結(jié)果及對比分析

分別計算27個模型的sig，R2和RMSE，圖4為3種建模方法所建模型的R2和RMSE對比結(jié)果，具體數(shù)值如表3所示。

(a)R2(b)RMSE

圖43種模型的R2和RMSE對比

Fig.4Comparisonofthreemodels’R2andRMSE

表327種估測模型檢驗樣本集結(jié)果統(tǒng)計

Tab.327hyperspectralestimationmodels’resultsoftestsamples

參數(shù)方法SQREEXPLOGESLSDE1LGD1EXD1均值R2BP0.8110.8070.8850.8120.8280.9330.8520.8390.8410.845MLR0.7780.8040.8350.7310.8590.8790.7920.8450.7860.812LSR0.8150.7970.8710.8260.8330.8870.8610.8470.8020.839RMSEBP0.1140.1170.0810.1100.1020.0690.0920.1070.1050.099MLR0.1220.1210.0990.1310.0930.0920.1180.1210.1230.113LSR0.1190.1250.0970.1140.1150.0850.0840.0920.1140.105sigBP0.0150.0190.0010.0110.0060.00 0.0030.0130.0070.008MLR0.0270.0310.0330.0420.0020.0010.0130.0060.0350.021LSR0.0190.0280.0010.0210.0050.0010.0020.0070.0160.011

由圖4(a)可以發(fā)現(xiàn)，BP與LSR模型的R2均較高，除幾種情況(RE和LGD1變換對2種建模方法效果相近； ES變換結(jié)合MLR模型效果優(yōu)于BP和LSR模型)外，均明顯優(yōu)于MLR模型，說明建模因子與土壤有機質(zhì)含量之間不是簡單的線性關(guān)系，MLR模型可靠性不高。同時，比較BP模型和LSR模型，發(fā)現(xiàn)對應(yīng)多種數(shù)據(jù)變換方法，二者所建模型的R2差別不大，但是LS變換對BP模型更敏感，所建模型的R2明顯優(yōu)于其他模型。由圖4(b)可以看出，BP模型的RMSE較小，MLR模型的RMSE較大，說明BP模型的估測精度高于MLR模型。綜合發(fā)現(xiàn)，BP模型結(jié)合LS變換所建模型的R2最大，RMSE最小，因此本實驗確定的最佳估測模型為BP+LS模型。表3中數(shù)據(jù)也顯示，所有模型的sig值均小于0.05，說明所建估測模型均是顯著的，結(jié)果可靠，所做比較具有意義。比較各模型R2和RMSE的均值結(jié)果可以發(fā)現(xiàn)，從整體而言，BP模型和LSR模型均優(yōu)于MLR模型，且BP模型更優(yōu)于LSR模型。

圖5為BP+LS模型實測值與估測值的對比結(jié)果，橫坐標(biāo)表示有機質(zhì)含量的實測值，縱坐標(biāo)為BP+LS模型得到的有機質(zhì)含量的估測值，由表3中數(shù)據(jù)可知，BP+LS模型的R2可達0.933，檢驗樣本的RMSE為0.069。

圖5 估測值與實測值對比結(jié)果Fig.5 Comparison of predicted and measured values

圖5中虛線為估測點的趨勢線，直線為y=x指標(biāo)線，在直線下方的點表明有機質(zhì)含量的估測值低于實測值，而直線上方的點表明有機質(zhì)含量的估測值高于實測值。從圖6可以看出，當(dāng)有機質(zhì)含量<0.6%時，估測值普遍高于實測值； 當(dāng)有機質(zhì)含量在0.6%～0.8%之間時，估測值大多低于實測值； 當(dāng)有機質(zhì)含量>0.8%時，估測值高于實測值。當(dāng)有機質(zhì)含量偏小或者偏大時，BP+LS模型得到的有機質(zhì)含量的估測值多高于實測值，但總體而言，差別不大，模型精度較高。LS變換有效提高了建模因子與有機質(zhì)含量的相關(guān)性，同時，BP神經(jīng)網(wǎng)絡(luò)具有較強的非線性映射能力，能較好地模擬土壤有機質(zhì)含量和地面實測高光譜反射率之間的關(guān)系，此外，BP神經(jīng)網(wǎng)絡(luò)在建模時，能通過學(xué)習(xí)自動實現(xiàn)輸出數(shù)據(jù)的合理規(guī)則，因此BP+LS估測模型具有較高的精度和可靠性。

3 結(jié)論

1)數(shù)據(jù)質(zhì)量對建模效果影響顯著。通過比較原始光譜數(shù)據(jù)和小波去噪后光譜反射率數(shù)據(jù)與有機質(zhì)含量的相關(guān)系數(shù)，發(fā)現(xiàn)經(jīng)小波分析處理后，在590～770 nm之間，光譜反射率與有機質(zhì)含量的相關(guān)性顯著提高。研究提高數(shù)據(jù)質(zhì)量及數(shù)據(jù)與目標(biāo)因子相關(guān)性的理論和方法具有重要意義。

2)比較不同數(shù)據(jù)變換的結(jié)果，發(fā)現(xiàn)數(shù)據(jù)變換可以提高數(shù)據(jù)和土壤有機質(zhì)含量的相關(guān)性，但是建模因子對不同數(shù)據(jù)變換方法響應(yīng)程度不同，不恰當(dāng)?shù)淖儞Q可能使參數(shù)的相關(guān)系數(shù)不增反降，綜合比較，對數(shù)的平方變換是一種有效的變換方法，對3種建模方法效果均較好，在同類模型中精度和可靠性最高，因為對數(shù)的平方運算能較好地消除數(shù)據(jù)異方差問題，縮小差距。

3)針對建模因子與有機質(zhì)含量相關(guān)性不高的問題，改變了以往對光譜反射率數(shù)據(jù)進行變換的思路，提出了對提取的建模因子直接進行數(shù)據(jù)變換，然后利用變換后的數(shù)據(jù)建立高光譜土壤有機質(zhì)含量估測模型的方法。本文利用土壤地面實測高光譜數(shù)據(jù)(350～2 500 nm和有機質(zhì)含量數(shù)據(jù)，通過小波分析和包絡(luò)線去除法提取特征參數(shù)，選取與有機質(zhì)含量相關(guān)系數(shù)較高的參數(shù)作為建模因子，經(jīng)過9種數(shù)據(jù)變換方法處理后，結(jié)合BP神經(jīng)網(wǎng)絡(luò)法、多元線性回歸法和最小二乘回歸法共建立27個土壤有機質(zhì)含量估測模型，對比各模型的評價指標(biāo)，最終確定BP神經(jīng)網(wǎng)絡(luò)+對數(shù)的平方變換模型為最佳估測模型。

本文探究了多種數(shù)據(jù)變換方法和建模方法的組合在高光譜土壤有機質(zhì)含量估測中的效果和精度，并最終確定了相對最優(yōu)的建模方法和估測模型，對該領(lǐng)域高光譜技術(shù)的應(yīng)用推廣具有一定的參考意義。但由于本文沒有顧及土壤水含量對有機質(zhì)含量估測的影響，仍存在一定的局限性，未來將繼續(xù)深入研究。