基于光譜指數的晉西黃土區土壤黏粒含量估測

孫瑞鵬，丁皓希，畢如田，鄧永鵬，朱洪芬

（山西農業大學 資源環境學院，山西 太谷 030801）

土壤黏粒含量是土壤最基本的物理性質之一， 對土壤的水肥氣熱等各種性質均有影響，土壤質地的實驗室測定方法主要有沉降法和激光粒度儀法，其中，激光粒度儀法出現時間短，相比沉降法所測定的黏粒含量偏低，測定過程比較費時且成本較高[1-5]。目前，已有學者采用高光譜技術預測土壤黏粒含量，常用光譜變換來提高土壤黏粒含量的預測能力。王德彩等[6]引入正交信號校正（Orthogonal signal correction，OSC）處理后的光譜對土壤黏粒含量進行預測，結果表明，預測結果OSC處理高于原始光譜和微分處理，且OSC處理能消除不相關因素的影響。張雅梅等[7]選用7種光譜變換形式對土壤黏粒含量進行預測，結果表明，原始光譜的對數為最佳的光譜變換形式，且對土壤黏粒含量的預測能力最佳。雖然使用預處理光譜能夠預測土壤黏粒含量，但這樣只是通過多波段信息來建立預測模型，預測所用光譜信息比較單一，并未結合較多的光譜特征信息。

在基于高光譜的預測模型方面，近年來較多采用非線性模型用于土壤屬性的預測研究。張娜等[8]采用一元線性回歸、逐步多元回歸和神經網絡預測土壤砂粒和粉粒含量，結果表明，神經網絡的預測效果更好。相比線性模型，采用非線性模型會提高土壤質地的預測精度，由于非線性模型的容錯性強且有較強的魯棒性，對預測模型起到了顯著的優化效果。另外，還有不少研究者通過構建光譜指數以提高其預測能力，并篩選出各種土壤光譜特征指數解釋其對光譜曲線的影響。例如，張娟娟等[9]構建差值光譜指數DI（CR1883，CR2065）較好地估測了土壤有機質、全氮和速效氮3種土壤養分。趙明松等[10]構建差值指數、比值指數和歸一化指數，結合弓曲差較好地實現了土壤有機質的預測。尼加提·卡斯木等[11]對光譜指數進一步優化處理，結果表明，土壤光譜指數能更好地預測土壤有機質。還有相關學者采用不同光譜指數對土壤鹽分、土壤電導率等土壤屬性進行預測研究得出，模型預測良好[12-13]。以上研究表明，土壤光譜指數可以較好地預測土壤的基本屬性。

黃土高原土壤顆粒組成主要以粉粒和砂粒為主，二者含量占顆?？偭康?0%以上[14]，土壤黏粒含量較低，抗蝕性較弱[15-16]。而晉西黃土區是水土侵蝕較嚴重的區域，土壤黏粒含量是評價水土侵蝕的重要因子，較低的土壤黏粒含量為其光譜預測帶來了挑戰。

本試驗以位于晉西黃土區的臨汾市大寧縣為研究區，基于土壤高光譜預處理數據構建光譜指數模型，通過相關性分析篩選出最優土壤黏粒光譜指數，并采用多種預測模型構建土壤黏粒含量預測的最優組合形式，旨在為晉西黃土區土壤黏粒含量的快速預測提供方法。

1 材料和方法

1.1 研究區概況

山西省臨汾市大寧縣位于黃河中游，晉西呂梁山南端，位置36°16′～36°36′N，110°28′～111°01′E。地形南北高、中間低，東部高、西部低，海拔最低481 m，最高1 740 m，屬暖溫帶大陸性氣候，年均氣溫11.1℃，年均降水量536.9 mm。根據國家土壤信息服務平臺（http://www.soilinfo.cn/map/）中國1∶400萬土壤類型圖可知，區域內土壤類型主要有黃綿土和褐土兩大類。

1.2 土壤樣品采集與質地測定

土樣采集時間為2020年11月，采樣深度為20 cm，所采土樣質量約1 kg，記錄采樣點的位置信息，共采集土樣192個（圖1）。土樣帶回實驗室內自然風干后，將土樣分成2份，一份原土過2 mm篩用于土壤質地測定；另一份研磨后過2 mm篩用于土壤高光譜測定。

土壤粒徑數據使用激光粒度儀Mastersizer 3000進行測定，測定粒徑范圍為0.02～2 000μm，測定前采用H2O2-HCl-（NaPO3）6法對土樣進行預處理：首先將2～5 g小于2 mm的土樣放入高型燒杯中，緩慢少量逐次加入30%的H2O2，直至沒有氣泡產生；然后多次滴入10%的HCl，直至碳酸鹽完全去除；最后多次使用蒸餾水洗至中性后，移除上清液并加入濃度為0.1 mol/L的（NaPO3）6分散劑，將樣品放到超聲波振蕩器中10 min，待樣品分散后用于測定[4，17-19]。根據美國農部制土壤質地分類法將測定所得土壤粒徑劃分為：黏粒（＜0.002 mm）、粉粒（0.002～0.050 mm）、砂粒（0.05～2.00 mm）。本研究選用土壤黏粒含量進行分析。

1.3 土壤光譜測定

土壤光譜使用ASD Field Spec 4 Std-Res地物光譜儀進行測定，波長范圍為300～2 500 nm，在黑暗的實驗室環境下進行。將土樣裝入深2 cm的培養皿中并壓實，光纖探頭采用10°鏡頭，保持與土樣垂直距離15 cm，光源天頂角為30°，距離土樣50 cm；每個土樣旋轉重復測定5次，以這5條光譜曲線的平均值作為土樣的光譜數據。由于數據在350～400 nm和2 450～2 500 nm這2段光譜波段的信噪比大，在計算分析中剔除，然后對400～2 450 nm光譜數據進行Savitzky-Golay濾波9點平滑處理后并按5 nm間隔進行重新采樣，共得到411個波段作為光譜分析數據。

1.4 光譜變換與光譜指數構建

首先將原光譜反射率（Reflectance，R）進行3種變換：倒數變換（Inverse of R，IR）、倒數的對數變換（Logarithm of 1/R，LGIR）、倒數的一階微分變換（First derivative of 1/R，FDIR）；然后將以上4種形式的光譜曲線通過光譜指數兩兩組合，主要有差值光譜指數（Difference spectral index，DSI）、比值光譜指數（Ratio spectral index，RSI）和歸一化光譜指數（Normalized difference spectral index，NDSI）；最后計算土壤黏粒含量，并將其與各個光譜指數進行相關性分析，在相關性等勢圖中用區域極值的方法，選擇區域內相關性較高的光譜指數建立模型。

式中，Rm和Rn分別代表m和n波段的土壤光譜反射率（m＞n）。

1.5 建立模型與精度評價

采用土壤黏粒含量和選取的土壤光譜指數建立模型，模型包括多元線性回歸（Multivariable linear regression，MLR）、偏最小二乘回歸（Partial least squares regression，PLSR）、反向傳播神經網絡（Back propagation neural networks，BPNN）；對土壤黏粒含量依次排序，按照3∶1將全部樣本分為建模集和驗證集，分別為144個和48個樣本；采用決定系數（Coefficients of determination，R2）、均方根誤差（Root mean squares error，RMSE）、相對分析誤差（Relative prediction deviation，RPD）3個參數對模型精度進行評價。

式中，yi為樣本i實測值，y^i為樣本i預測值，yˉ為樣本實測平均值，n為樣本數，SD為樣本實測值的標準差。

R2值越大、RMSE值越小，表示模型精度越高。本研究將RPD值分為6類：RPD＜1.0表示模型預測非常差，不推薦使用；RPD在1.0～1.4表示模型預測比較差，只能用于區分高值和低值；RPD在1.4～1.8表示模型可用于預測；RPD在1.8～2.0表示模型預測良好；RPD在2.0～2.5表示定量模型預測非常好；RPD＞2.5表示模型預測極好[20-21]。

1.6 數據分析

光譜數據預處理使用View Spec Pro 6.2軟件完成；光譜變換、光譜指數、相關性計算和模型計算均使用Matlab 2019b軟件完成。

2 結果與分析

2.1 不同土壤黏粒含量的光譜特征分析

根據美國農部制土壤質地分類法[22]，192個土壤樣本中有188個為粉壤土，3個為粉土，1個為砂壤土，其中，土壤黏粒含量最小值為4.03%，最大值為17.83%，平均為8.00%；按照土壤黏粒含量從小到 大 排 序 為1～17、18～64、65～120、121～165、166～184、185～192個 土 樣 的 平 均 值 分 別 是5.00%、6.32%、7.50%、8.89%、10.85%、16.02%，共6組，計算每一組土樣的光譜平均值如圖2-A所示。土壤黏粒含量對土壤光譜的影響主要體現在近紅外波段；隨著土壤黏粒含量的增加，1 450 nm和1 950 nm處的水分吸收明顯，即光譜吸收谷越深（圖2-B、C），土壤黏粒含量越高，光譜曲線在近紅外波段反射率越低，且光譜差異較大，反之，光譜曲線在可見光波段差異較小。AL-ABBAS等[23]在對土壤質地與光譜的研究中指出，土壤黏粒含量與其光譜反射率呈負相關，土壤黏粒含量越高反射率越低，與本研究結果一致。

2.2 土壤黏粒含量與光譜指數的關系

通過分析土壤黏粒含量與原始光譜（R）、倒數（IR）、倒數的對數（LGIR）及倒數的一階微分（FDIR）的相關性可知（圖3），相關系數較高值分別為-0.28、0.27、0.28、-0.61，經IR和LGIR變換后的相關系數并未提升，而經FDIR變換后其相關性明顯提升，相關性較高的波段主要位于1 450、1 950、2 200 nm光譜吸收谷兩側。

圖4為光譜變換對應的差值（DSI）、比值（RSI）及歸一化（NDSI）光譜指數與土壤黏粒含量的相關關系情況。從圖4可以看出，在同一種光譜變換形式下，DSI、RSI和NDSI這3種指數的正負相關性區域相近，其中，DSI相關性圖與RSI和NDSI存在較大差異，而RSI和NDSI相關性區域基本保持一致；在同一種光譜指數下，IR和LGIR這2種變換正負相關性區域比較吻合，且二者與原始光譜的正負相關性區域正好相反。總體來看，相關性區域相對集中，相關性最高的區域分布在波段1 450、1 950、2 200 nm吸收谷的兩側，波段組合為（2 245，1 515）、（2 425，1 155）、（2 440，1 155）、（1 900，1 460）、（2 240，1 515），相關系數絕對值在0.66～0.70；反之，FDIR變換對應的相關性區域比較分散，DSI、RSI、NDSI對應的顯著相關性最高波段組合分別為（1 930，1 420）、（2 420，1 570）、（1 940，400），相關系數分別為-0.67、-0.70、0.69；相較于單波段，波段兩兩組合與土壤黏粒含量的相關性均有提升。從土壤黏粒含量與光譜指數的相關性來看，光譜指數RSI最優，NDSI和DSI表現次之（表1）。

由于相關系數區域比較集中，為減少構建模型時自變量之間共線性的影響，對相關系數圖進行求極值，然后在每種光譜指數中篩選出相關系數極值絕對值最大的5個土壤光譜指數?；谝陨献顑炌寥拦庾V指數，采用多元線性回歸方法構建土壤黏粒含量的預測模型，建模結果如表1所示，12種形式的建模集R2在0.47～0.57，RMSE在1.42～1.59，對應大小與每一類型相關系數絕對值大小有關，波段組合與土壤黏粒含量相關性越高，建模精度越高；驗證集R2在0.45～0.61，RMSE在1.42～1.71，RPD在1.30～1.57，除FDIR-DSI和FDIR-NDSI之外，其余10種形式的RPD均大于1.4，能夠粗略地估測土壤黏粒含量；其中預測精度最高的是FDIR-RSI，建模集的R2和RMSE分別為0.57、1.42，驗證集的R2、RMSE和RPD分別為0.61、1.42、1.57，其 次 分 別 為R-DSI、R-NDSI、IRNDSI、LGIR-DSI。

表1 各光譜指數與土壤黏粒含量的關系Tab.1 Relationship between each spectral index and soil clay contents

2.3 基于PLSR和BPNN的土壤黏粒含量估測

與MLR采取同樣的光譜指數選取方法，通過極值結果選取相關系數絕對值最大的47個光譜指數，用于PLSR和BPNN預測，模型結果如表2所示，PLSR模型建模集R2和RMSE分別在0.60～0.67、1.25～1.38，驗證集R2、RMSE和RPD分別在0.50～0.57、1.46～1.56、1.42～1.52，其中，RPD均大于1.4，預測精度高的是FDIR-RSI，建模集的R2和RMSE分別為0.63、1.32，驗證集的R2、RMSE和RPD分別為0.57、1.54、1.52，預測精度其次依次為R-DSI和LGIR-RSI；BPNN模型建模集R2和RMSE分別在0.64～0.70、1.21～1.32，驗證集R2、RMSE和RPD分別在0.55～0.74、1.13～1.49、1.50～1.96，同樣RPD均大于1.4，預測精度高的是LGIRNDSI，建模集的R2和RMSE分別為0.64、1.32，驗證集的R2、RMSE和RPD分別為0.74、1.13、1.96，預測精度其次為LGIR-RSI、LGIR-DSI和FDIRRSI。總體來看，在LGIR變換下3種光譜指數構建的BPNN模型精度都比較好，整體比較DSI、RSI和NDSI對應的模型精度，RSI和NDSI優于DSI。比較PLSR、BPNN模型建模集和驗證集的R2和RMSE可知，BPNN表現最優，其次為PLSR。

表2 基于光譜指數的土壤黏粒含量估測結果Tab.2 Estimation of soil clay contents based on spectral indexes

圖5是模型預測效果，對比實測值與預測值分布，存在低值高估和高值低估的現象；BPNN預測結果更接近1∶1線，預測結果有所改善，非線性模型BPNN用于預測土壤黏粒含量比傳統線性模型表現更好。

3 結論與討論

已有研究采用光譜指數預測其土壤屬性，模型預測使用光譜范圍主要為近紅外光譜區域。同樣，本研究結果表明，土壤黏粒含量與光譜指數顯著相關性的區域主要分布在近紅外區域，主要為（2 245，1 515）、（2 425，1 155）、（2 440，1 155）、（1 900，1 460）、（2 240，1 515）等相鄰區域，這與已有光譜指數預測土壤屬性結果一致。郭鵬等[13]通過構建光譜指數預測土壤鹽分，發現相關系數高的光譜范圍為1 430～1 862、1 934～2 150 nm，篩選出敏感亮度指數（1 750，1 620）,并采用隨機森林（Random forest，RF）建立預測模型，且模型精度較高。張娟娟等[9]構建差值光譜指數DI（CR1883，CR2065）用于估測土壤有機質、全氮和速效氮3種土壤養分，結果發現均有較好的效果。

在光譜變量選取方面，本研究采用光譜指數結構固定、光譜波段隨機組合的方式構建光譜指數，并使用與土壤黏粒含量的相關性篩選局部區域的最優光譜指數，以提升土壤黏粒含量的預測精度。已有研究常采用相關性、主成分變換等方法減少輸入變量，例如，楊莎等[24]使用相關性分析篩選敏感波段用于預測土壤有機質。白燕英等[25]通過土壤黏粒含量與光譜的相關性選擇波段用于建立預測模型，精度超過85%。喬天等[26]采用遺傳算法選擇變量并結合偏最小二乘回歸建立土壤質地模型，結果減少了光譜變量，也提高了預測精度。王德彩等[27]采用光譜主成分為變量使用BP神經網絡預測土壤黏粒含量，實現了高效預測。而本研究使用光譜指數一方面避免了篩選出相鄰光譜存在的光譜互相關，提高了預測模型的穩健性；另一方面光譜指數能夠更多選取到不同的光譜波段，涉及的光譜信息更多。

在土壤光譜預測方面，較多使用PLSR模型預測土壤屬性。本研究同時使用PLSR和BPNN對土壤黏粒含量進行預測，結果表明，BPNN預測能力更強更穩定。同樣，郭云開等[28]使用BPNN預測土壤Cu含量，結果表明，相比單元線性回歸模型，BPNN具有良好的擬合優度和預測能力。雖然本研究使用光譜指數預測了晉西黃土區的土壤黏粒含量，可以為遙感估測土壤黏粒含量提供技術支持，但相比較大尺度范圍遙感估測，研究區域還較小，土壤類型較單一，導致光譜指數預測土壤黏粒含量模型適用性較差，所以該方法用于不同區域還有待進一步探究。

本試驗基于光譜指數的晉西黃土區的土壤黏粒含量估測研究結果顯示，土壤黏粒含量對光譜曲線的影響主要為近紅外波段，原始光譜與土壤黏粒含量呈負相關，土壤黏粒含量越高，土壤光譜反射率越低；通過構建光譜指數可以提高與土壤黏粒含量之間的相關性，是提高預測土壤黏粒含量的方法之一；基于LGIR-NDSI-BPNN的預測效果最好，驗證集R2和RPD分別為0.74和1.96，能夠有效地預測晉西黃土區的土壤黏粒含量。為了提高光譜指數預測土壤黏粒含量的精度，將針對不同土壤類型、多種光譜指數探究最佳的光譜預測指數和方法。