基于連續小波變換的土壤有機質含量高光譜反演

陳昊宇, 楊光, 韓雪瑩, 劉昕, 劉峰, 王寧

(內蒙古農業大學沙漠治理學院， 內蒙古自治區風沙物理與防沙治沙工程重點實驗室， 呼和浩特 010010)

精準農業作為目前農業發展的主要方向，是一種基于信息和知識管理的現代化生產系統，主要是通過3S(GPS、GIS和RS)技術與現代農業相結合，最大限度地提高農業生產力。所以快速、無損、精確地獲取土壤中水分、養分的空間分布成為了實現精準農業的關鍵環節，近年來，光譜分析在土壤化學分析領域得到了迅猛發展，為實現土壤養分的快速診斷提供了新思路[1]。有機質是土壤養分供應能力和肥力的重要指標之一，在全球碳循環中發揮著重要作用。因此，快速準確地估測土壤有機質含量對于發展精準農業具有重要意義[2]。

傳統的土壤有機質測定方法雖然精度比較高，但周期較長、成本較高，只能達到瞬測量，很難進行長時間大面積測量。高光譜遙感具有波段多、波段窄、信息豐富和實時高效等特點，為快速測量土壤有機質含量提供了一種新的方法和手段[3]。

國內外已經有大量研究表明，通過對光譜數據進行不同的數學變換(主要通過對光譜進行倒數、對數、微分、平方根、吸收峰深度、包絡線去除等方法)可以有效提高光譜數據與土壤有機質含量之間的相關系數，有效篩選出光譜信息中的敏感波段[4]?，F在各學者主要將研究重心放到了模型建立上[5]，普遍運用的線性模型有多元逐步回歸與偏最小二乘回歸[6]；常見的非線性模型包括BP神經網絡[7]、支持向量機[8]、決策樹[9]等，而且隨著非線性模型算法的逐步改良與完善，在土壤有機質含量估算中已經成為不可取代的一部分。隨著小波算法的改進與發展，最初僅運用于植物葉綠素、冠層成分含量預測中[10-11]，目前已成為土壤養分預測的熱點問題[12-13]，連續小波變換是目前被廣泛應用的一種方法。王祥浩[14]選擇土地裸露地區為樣區，利用神經網絡算法對光譜連續小波變換、一階導數、對光譜的平均值處理、光譜背景及深度4種方法建模，模型結果表明，小波變換方法得到的神經網絡模型精度最高；包青嶺等[15]選擇渭干河-庫車河三角洲具有代表性的干旱區綠洲為研究區，對光譜進行8層分解，結果表明小波變換不同分解層，從低頻到高頻范圍內與土壤有機質含量的相關性呈現先減后增的趨勢，結合隨機森嶺模型可以對干旱區土壤有機質含量進行有效的估算；王延倉等[16]以北京東部區潮土為例，對不同梯度重采樣的光譜進行連續小波變換后，利用偏最小二乘法建立模型，結果表明連續小波分析算法可深入挖掘土壤光譜內的有益信息，提升對有機質含量的估測能力，與土壤高光譜反射率相比，經連續小波技術處理后，模型精度得到了有效的提升；葉紅云等[17]同樣針對干旱區土壤，通過對兩種常用光譜變換R′、Ln(1/R)進行連續小波變換建立偏最小二乘模型，結果表明連續小波變換不會因人類干擾程度的提高而使模型精度大幅度降低，更加適用于干旱區有機質含量的預測；林鵬達等[18]通過解決黑土有機質高光譜野外反演的困難，同樣證明了連續小波變換可有效提升模型精度。小波技術在土壤有機質高光譜反演研究中逐漸趨于成熟，但目前學者的研究多數都在同一土壤類型下或同一區域內，對于不同土壤類型及土地利用下土壤有機質高光譜反演是否存在影響的研究目前并不多。本文研究區內土壤類型主要包括3類：沙壤土、栗鈣土、鹽堿土，且部分區域土壤鹽漬化程度嚴重，導致土壤養分空間分布上存在較大差異，取樣表層土地利用類型主要包括：耕地、林地、草地、鹽漬地、荒地。

通過對原始光譜(R)、原始光譜倒數(1/R)、原始光譜對數(LnR)以及原始光譜一階微分(R′)4種不同情況進行連續小波變換，利用BP神經網絡以及支持向量機2種模型，探究了不同土壤類型與不同土地利用類型下是否會對土壤有機質高光譜反演模型產生影響，小波變換前后土壤有機質反演模型的精度，旨為區域土壤有機質含量監測及實現精準農業提供理論與技術支持。

1 材料與方法

1.1 研究區概況

托克托縣隸屬于內蒙古自治區呼和浩特市，位于自治區中部、大青山南麓、黃河上中游分界處北岸的土默川平原上(圖1)。地理坐標東經111°2′30″—111°32′21″、北緯40°5′55″—40°35′15″，總面積1 409.67 km2，平均海拔1 117 m，屬于溫帶大陸性干旱氣候，年均氣溫7.3 ℃，年均降雨362 mm。托克托縣耕地總面積達400 km2，其中古城鎮、新營子鎮和五申鎮的耕地較多，占全縣耕地面積的60%以上[19]，主要作物包括小麥、玉米、莜麥。工農業及生產生活用水主要來源于大黑河和黃河水資源，整個地形以大黑河為軸，呈現由丘陵向平原過渡的趨勢，地勢為東南高、西北和西南低。東南向西北土壤類型依次為栗鈣土、砂壤石灰性沖積土、鹽漬化石灰性沖積土[20]，土壤類型的不同導致土壤養分存在差異性分布。植被類型從西向東依次為草甸草原、干草原和退化灌叢草原分布。以Landsat8OLI影像為基礎數據源，運用人工目視解譯與BP神經網絡分類法得到托克托縣2019年7月份土地利用數據，其中耕地面積最大為730.12 km2，占51.79%；林草地338.7 km2，占24.02%；鹽堿地141.1 km2，占10.00%。詳細土地利用空間分布見圖1。

圖1 土樣采集點及土地利用空間分布

1.2 研究方法

1.2.1土樣采集與處理土壤樣本點均勻地分布在托克托縣境內，采集方法為五點采樣法，采集深度為0—20 cm，共采集120個點。采集的土樣置于通風干燥室內進行自然風干、研磨，過10目篩，進行土壤光譜測定；過100目篩，采用重鉻酸鉀外加熱法進行土壤有機質含量測定。

1.2.2光譜測量及光譜處理土壤光譜于暗室內測量，采用SVC HR-1024(北京東方佳氣科技有限公司)便攜式光譜儀，光譜范圍在350～2 500 nm。在350～1 000 nm波段之間光譜分辨率≤3.5 nm；在1 000～1 850 nm波段之間，光譜分辨率≤9.5 nm；在1 850～2 500 nm波段之間，光譜分辨率≤6.5 nm。光源采用與太陽光接近的50 W鹵素燈，將土壤樣品放入深2 cm、寬10 cm的黑色器皿內，用直尺將土壤表面刮平，探頭距離土樣10 cm，光源距離土壤表面30 cm，天頂角為15°。測量前用白板進行標定，每個土樣采集5條光譜作為該土樣的光譜數據。

由于受噪音與儀器暗電流的的影響，導致光譜數據混入噪音等信息，因此刪除350～399 nm和2 400～2 500 nm的波段，采用五點平滑法對光譜進行平滑處理，并將光譜重采樣至5 nm，同時對原始光譜(R)進行一階微分(R′)、倒數(1/R)、對數(LnR)等傳統數學變換。

1.2.3連續小波變換采用連續小波變換，并用Mexh小波母函數對原始光譜、原始光譜的倒數、對數、一階微分進行10層小波變換，生成一系列小波系數。

(1)

式中，a為伸縮因子，b為平移因子，λ為土壤高光譜數據的波段數。

(2)

式中，f(λ)為土壤光譜反射率，小波系數Wf(a，b)包含二維，分別為波長(350～2 500)與分解尺度(1，2，3…10)， 故小波系數行為尺度數，列為波長數的矩陣[16]。

1.2.4模型及精度驗證采用BP神經網絡與支持向量機模型(support vector machine，SVM)建立土壤有機質預測模型，支持向量機采用線性核函數，相對于徑向基函數(radial basis function, RBF)來說計算高效,不易過擬合。BP神經網絡的迭代次數設置為1 000，學習率0.01，訓練的均方根誤差(root mean square error，RMSE)小于0.001。

依據相關系數篩選的特征波段以及小波系數作為自變量，土壤有機質含量為因變量，分別建立模型，模型精度采用決定系數(R2)、均方根誤差(RMSE)、相對分析誤差(relative percent deviation，RPD)以及1∶1線共同評價。R2表征模型的穩定性，越接近于1模型越穩定，擬合程度越好。均方根誤差(RMSE)用來檢驗模型的預報能力，RMSE越小則表明模型的估測能力越好。RPD是樣本的標準差與RMSE的比值，RPD<1.4時，模型無法對樣品進行預測；1.4≤RPD<2時，模型效果一般，可以用來對樣品進行粗略評估；RPD≥2時，模型具有極好的預測能力。1∶1線表示實測值與預測值構成的點偏離y=x線的程度[21]。

2 結果與分析

2.1 土壤有機質含量統計分析

建模樣品集、不同土地利用方式、不同土壤類型下土壤有機質含量描述性統計見表1。本研采樣點內土地利用方式主要包括林地、草地、耕地、鹽漬地，土壤有機質在草地內均值含量最大(0.80%)，其次為林地(0.72%)、耕地(0.67%)、鹽漬地有機質含量最低(0.63%)；土壤有機質含量最大值位于耕地(1.28%)，最小值位于林地(0.19%)。采樣點內主要土壤類型為栗鈣土、沙壤土、鹽堿土，沙壤土有機質含量最高(0.77%)，其次為鹽堿土(0.68%)和栗鈣土(0.67%)，土壤有機質含量最大值位于沙壤土內(1.28%)，最小值位于鹽堿土內(0.19%)。

表1 土壤有機質含量描述性統計結果

2.2 土壤反射光譜特征

對R、1/R、LnR、R′進行小波變換，變換結果如圖2所示，R、1/R、LnR光譜曲線較為平滑，分解曲線隨波峰波谷變化.R′其光譜曲線并不規則存在較多波峰波谷，分解小波系數與前三者不同。R、1/R、LnR、R′分解后，小波系數均隨分解尺度的增加而增加，同時可以看出，由Mexh小波母函數進行的連續小波變換，對于光譜波峰與波谷有較高的敏感性，對于放大、挖掘光譜信息有著顯著的作用。

圖2 連續小波變換光譜特性

2.3 相關性分析

2.3.1不同導數變換光譜與土壤有機質含量相關性土壤有機質含量與光譜相關性曲線及敏感波段見圖3。R與土壤有機質含量呈負相關關系(相關系數r=-0.463)，主要集中于735～780 nm處波段；1/R與土壤有機質的相關性則與R相反，呈正相關關系(r=0.462)，集中于600～800 nm與1 800～2 200 nm處波段；LnR的相關性曲線圖與R相關性曲線類似，總體呈現負相關關系，相關系數(r=-0.465)，主要集中于745～795 nm處的波段；R′相關性在500 nm(r=-0.589)與1 400 nm(r=-0.411)處為負相關，在8 00 nm(r=0.408)與1 380 nm(r=0.412)處為正相關，相關系數曲線變換趨勢與前三者不同，呈無規律變化。

圖3 土壤光譜相關性曲線及敏感波段

2.3.2不同分解尺度小波系數與土壤有機質含量的相關性圖4為不同光譜變換方式經過連續小波變換后與土壤有機質含量的相關系數矩陣圖，其中紅色代表相關性高的區域，藍色代表相關性低的區域。R在800～1 000、1 400～1 600 nm處相關性明顯增加，在500、800、2 200 nm波段處相關系數達到最大值(r=0.667)；1/R在800～1 200 nm處相關系數達到最大值(r=0.552)，在2 400～2 500 nm處相系數達到0.4，受噪音和儀器本身的影響，此波段的相關系數不進行相關性參考；LnR在分解尺度1下相關性較低，在2～10尺度下，相關性出現最大值(r=0.664)；R′相關性主要集中在500～900、1 200～1 600、2 100～2 300 nm處。篩選的敏感波段與尺度如表2所示。有效的光譜信息主要存在于低分解尺度，隨分解尺度的增加呈遞減趨勢，相關性最大值較未處理前分別增加了0.204、0.09、0.199、0.252，對于挖掘潛在光譜信息有著重要意義。

圖4 土壤有機質與小波系數相關性

表2 篩選的敏感波段

2.4 土壤有機質高光譜模型建立

2.4.1BP神經網絡預測模型采用BP神經網絡構建反演模型，結果如表3所示。未進行連續小波變換處理的模型中，BP-R與BP-R′效果較好，R2分別為0.69和0.73，RPD為1.45與1.53，模型能粗略估算土壤有機含量，BP-LnR與BP-1/R樣本外預測能力較差，同時RPD未達到1.4以上，不能對土壤有機質未能進行有效預測；連續小波變換處理之后的模型，僅BP-CWT-1/R模型RPD未達到預測水平，其余3種模型R2與RPD較未處理前均有所增加，RMSE均減少，其中BP-CWT-LnR模型預測效果較好，RPD達到2.12可以有效地對土壤有機質進行預測。將BP-CWT處理的4個模型的實測值與預測值進行1∶1線分析。由圖5可知，除BP-CWT-1/R模型外，其余模型的實測值與預測值樣點基本分布在1∶1線附近，BP-CWT-LnR效果較為明顯，且估算精度高，可較好地進行土壤有機質含量的估算。

表3 土壤有機質BP神經網絡估測模型結果

圖5 BP-CWT模型土壤實測值與預測值對比

2.4.2支持向量機預測模型SVM構建反演模型，結果如表4所示。未經過連續小波處理的光譜特征波段未能較好地對土壤有機質進行預測反演，經過CWT后模型SVM-CWT-R與SVM-CWT-R′預測結果較之前有較大的提升，R2分別達到了0.50與0.56，二者RPD均達到1.4以上，可以粗略地對土壤有機質進行預測。同時根據圖6，SVM-CWT模型進行1：1線分析，二者實測值與預測值分布情況在4種模型下較好，雖然模型SVM-CWTLnR分布同樣較為集中，但其樣本外預測情況較差(RPD=1.38)，綜合考慮不對其進行土壤有機質預測。結合表3和表4的結果分析，連續小波變換能夠有效地提升模型精度與模型泛化能力，對于光譜信息挖掘有著重要意義，BP神經網絡與支持向量機對CWT-R與CWT-R′都能夠提升R2減少RMSE，可對土壤有機質做出較好的預測。雖然BP神經網絡與支持向量機在處理非線性回歸問題中有較強的能力，但本身模型中存在不穩定性，對模型的環境設置同樣要求較高，所以未能對所有數據集進行良好的預測。

表4 土壤有機質支持向量機估測模型結果

圖6 利用SVM-CWT模型土壤實測值與預測值的對比

3 討論

本研究采用連續小波變換對光譜進行處理，用BP神經網絡與支持向量機(SVM)兩種模型對土壤有機質含量進行反演預測。未經過連續小波變換前，R、1/R、LnR、R′與土壤有機質的相關系系數最大值分別為-0.463、0.462、-0.465、0.589，可以看出，R′與土壤有機質的相關系數最高，與吳倩等[22]、張新樂等[23]的研究結果相同；經過連續小波變換后，CWT-R、CWT-1/R、CWT-LnR、CWT-R′相關系數最大值分別為0.667、0.552、0.664、0.662，較之前分別增加了0.20、0.09、0.19、0.07。王延倉等[1]、于雷等[4]、葉紅云等[17]等同樣證明連續小波變換可有效提高與土壤有機質含量的相關系數。不同分解尺度對于光譜數據的深度挖掘有著重要意義，本研究只利用Mexh小波母函數進行處理，未對其他函數進行考慮，分解層數同樣是根據前人經驗所得[4，10]，小波技術的研究與發展仍然有很大的探索空間。

相對于兩種模型來看，未進行連續小波處理的支持向量機模型中，只有SVM-R′模型R2最高達到0.43，其余三者均未到達0.4。綜合多種模型評價方法，由于其RPD未達到1.4以上，無法對土壤有機質含量進行預測。經過連續小波處理后，各模型的R2有明顯提高，其中SVM-CWT-R與SVM-CWT-R′模型效果較好，R2分別提高了0.29、0.13，RPD達到1.62與1.53實現了對土壤有機質有效的預測，但預測結果較BP神經網絡較低。在BP神經網絡預測模型中，未進行連續小波變換前，BP-R與BP-R′預測效果較好，R2達到0.69與0.73，RPD為1.45與1.53；進行連續小波處理后，除SVM-CWT-1/R模型未到達預測效果，其余3種模型預測結果較之前均有明顯改善，可實現對土壤有機質較好的預測，其中BP-CWT-LnR預測模型效果最佳R2達到0.76，RPD達到2.12。根據1：1線分析圖也可看出，其實測值與預測值分布較為集中，于雷等[4]、葉紅云等[17]、林鵬達等[18]同樣通過連續小波變換有效提升了模型的精度與泛化能力。

針對土壤有機質高光譜反演研究中，姚聰[24]對耕層土壤通過BP神經網絡與支持向量機模型，反演精度R2分別為0.42與0.67；葉紅云等[17]采用連續小波變換對干旱區土壤有機質反演，模型精度R2=0.75、EMSE=0.71；謝文[25]在森林土壤有機質反演研究中，BP神經網絡模型R2=0.78、EMSE=0.77，支持向量機模型R2=0.87、EMSE=0.76。本研究對耕地、林草地、鹽堿地、栗鈣土、沙壤土、鹽漬土等不同土地利用類型與土壤類型進行綜合反演，最佳反演模型為BP-CWTLnR，R2=0.76、EMSE=0.15、RPD=2.12，與前人研究的結果基本相符，證明通過連續小波變換處理，不同土壤類型與土地利用類型未對土壤反演模型精度產生影響。所以采用連續小波變換進行光譜數據挖掘，采用BP-CWT-LnR神經網絡建立反演模型，可對不同土地利用于土壤類型條件下土壤有機質高光譜反演提供一定的理論支持與應用價值。