基于高光譜技術的退耕還林地年限判別

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

（山西農業大學資源環境學院，太原 030000）

0 引 言

黃河中游的大部分支流均處于黃土高原區，面積約為34.4×10km，是中國分布面積較大的生態脆弱區。為改善黃河中游的生態環境問題，中國于1999年試行，并于2000年正式開始實施退耕還林（草）工程。隨著退耕還林工程長久實施，研究不同年限退耕還林土壤的綜合性質對評估退耕還林工程的效果具有重要意義。由于部分地區無法提供有效的退耕資料，因此快速、準確地識別退耕還林地的年限對評估退耕還林工程的生態效益具有重要價值。

近年來，光譜分析技術被廣泛應用到中藥材種植年限識別、儲藏年限識別，野生食用菌儲藏年限識別、產地識別，農作物種類區分、茶樹種類識別、土壤類型分類等不同領域內。陳澤炎等對目前中藥材生長年限鑒定方法進行了總結，發現光譜分析技術可以高效、快速地識別不同種植年限的中藥材。卜海博等采用近紅外光譜法對不同生長年限的林下山參和園參進行種類和年限識別，結果表明近紅外光譜法準確可靠、快速無損。楊天偉等采用傅里葉變換紅外光譜技術，準確區分了不同產地、不同年份美味牛肝菌樣品，為野生食用菌的鑒別分類提供了技術參考。虞佳維等以黃淮海地區玉米、小麥和楊樹三種主要植被為研究對象，采用高光譜技術進行區分，結果表明光譜技術可有效用于不同植被的區分。翁海勇等采用近紅外光譜分析技術實現了茶樹葉片表沒食子兒茶素沒食子酸酯（Epigallocatechin Gallate ，EGCG）含量的快速、準確識別，從而為高EGCG含量茶樹的育種提供技術支持。史舟等基于中國西藏、新疆、黑龍江、海南等地的16種土類的1 581個土壤樣本，以光譜反射率一階微分主成分數據為自變量進行土壤光譜分類，并根據分類結果建立了中國不同土壤類型的有機質預測模型。趙小敏等以江西紅壤地區的7種紅壤亞類土壤為研究對象，對不同亞類及土屬的土壤高光譜特征進行分析，選取特征變量進行Fastclus聚類分析，結果表明，土壤亞類分類的準確度為86.23%，土屬的準確率僅為66.37%。目前，尚未發現開展退耕還林地不同年限識別的研究，應用光譜分析技術可以為不同年限退耕還林地的分類提供技術參考。

本研究以黃河中游大寧縣不同年限退耕還林地土壤為研究對象，以不同退耕年限刺槐林地土壤光譜曲線和多種土壤理化數據為基礎，分析退耕年限對土壤光譜及其綜合屬性的影響。利用倒數的對數（Reciprocal of Logarithm，RL）、一階微分（First Order Differential，FD）等光譜處理方法，結合主成分分析（Principal Component Analysis，PCA）、K 均值聚類（K-means Clustering Algorithm，K-means）、支持向量機（Support Vector Machine，SVM）和線性判別分析（Liner Discriminant Analysis，LDA）等方法探索不同退耕年限土壤的分類結果，實現基于土壤光譜的退耕年限快速區分，以期為快速準確識別退耕年限提供技術參考。

1 材料與方法

1.1 研究區概況

本研究以山西省大寧縣境內退耕還林地為研究區（110°27′55″～111°0′40″E， 36°16′40″～36°36′25″N）。大寧縣地處黃河中游，地勢南北高、中間低，屬典型的黃土殘垣溝壑區，是黃河泥沙的主要輸入區之一。2000年大寧縣被定為黃河中上游退耕還林（草）生態工程試點縣之一，配置有刺槐、檸條等樹種。查詢國家土壤信息服務平臺(http://www.soilinfo.cn/map/)中國1:400萬土壤類型圖可知，研究區內主要土壤類型為黃綿土和褐土。

1.2 土壤樣品采集與處理

2020年11月在研究區內分別選擇退耕3 a（RF）、8 a（RF）、14 a（RF）和18 a（RF）的刺槐林地各32處，共128個樣點（圖1a）。由于黃綿土土質疏松，易受侵蝕，水土流失嚴重，是實施退耕還林工程的主要土壤類型，因此樣點土壤類型均為黃綿土。土壤樣品采集選用“S”形布點法，采集0～20 cm表層土壤，每個采樣點取5個點的土壤混合后作為該采樣點土樣。將土樣分為兩份，一份土樣在室內自然風干、研磨并過2 mm篩，用于光譜測試，土壤有機碳SOC含量測試；一份原狀土樣過2 mm篩用于土壤質地測定。同時用環刀進行取樣，測定飽和導水率、土壤含水率。其中SOC含量測定采用重鉻酸鉀氧化—外加熱法測定，土壤含水率采用經典烘干法測定，飽和導水率利用恒定水頭法測定。土壤質地數據測定采用激光粒度儀Mastersizer 3000進行測定，測定粒徑范圍為0.02～2000m，測定所得土壤粒徑數據根據美國農部制粒級分級劃分為：黏粒（＜0.002 mm），粉粒（0.002～0.05 mm），砂粒（0.05～2 mm）。

1.3 光譜測試及數據處理

采用美國ASD（Analytica Spectra devices，Inc）公司的FiledSpec 4便攜式光譜儀測定土壤反射光譜曲線，其波譜范圍為350～2 500 nm。土壤光譜測量在暗室中進行，每次測量開始前將光譜儀預熱2 h。將土壤樣本放于直徑為8 cm，深度為2 cm的培養皿中，用直尺輕輕刮平表面，減小土壤表面粗糙度。光源由一個50 W的鹵素燈提供，光纖探頭視場角為10°，距離土樣表面15 cm。每次測試前采用標準白板校正，每個土樣測定10條光譜曲線，取平均值作為該土樣的實際光譜曲線。

圖1 研究區概況及土壤光譜特征參數Fig.1 General situation of the study area and soil spectral characteristic parameters

去除噪聲較為強烈的350～399、2 451～2 500 nm波段。同時，為避免數據冗余，對原始光譜反射率數據（Reflectance，）做10 nm重采樣，并對其進行Savitzky-Golay平滑（Savitzky-Golay smooth，SG），對平滑后的光譜數據分別做一階微分、倒數的對數、去包絡線（Continuum Removal，CR）和光譜特征參數（Spectral Characteristic Parameter，SCP）處理，最后對、FD、RL等數據進行主成分分析。其中，CR又稱為連續統去除法，可以將光譜反射率歸一化到0～1之間，有效突出光譜特征。SCP包括：吸收谷（）、面積（）、深度（DP）、吸收波段波長位置（）、斜率（）、寬度（）、對稱度（）等（圖1b），對稱度=/。相關數據處理在The Unscrambler X 10.4、ENVI 5.3、MATLAB R2019b以及SPSS 26.0中完成。

1.4 分類指標選取

分析不同年限退耕還林土壤光譜特征并進行試驗，最終以原始反射率主成分（R-PCA）、倒數的對數主成分（RL-PCA）、一階微分主成分（FD-PCA）、去包絡線主成分（CR-PCA）、SCP作為分類指標，其中SCP包含：第1、2、3個吸收谷的深度（DP、DP、DP），第1、2、3個吸收谷的對稱度（、、），第一個吸收谷的面積及前兩個吸收谷面積之和（、+），495～595、895～995、1 450～1 600 nm這3個波段去包絡線的斜率（、、）。

1.5 模型建立與精度評價

本研究采用K-means、SVM、LDA三種方法構建不同年限退耕還林地土壤光譜分類模型，對比分析三種模型的分類精度。K-means優點是處理大樣本時精度較高，本研究先對指標標準化以消除指標量級對結果的影響。SVM對解決中小樣本、非線性和高維數據具有獨特的優勢，本研究選用徑向基核函數，懲罰參數（）和Gamma參數通過交叉驗證方法進行選擇，為3，Gamma參數為0.33；LDA常用于降維和分類，核心思想是投影后類內方差最小，類間方差最大，本文LDA變換的時候主要是對自變量數據即光譜數據進行了變換，重新投影到一個超平面上，使得樣本有更大的類間距離，更小的類內距離。分類模型構建分別在SPSS 26.0、MATLAB R2019b和The Unscrambler X 10.4中完成。

將不同退耕年限土壤樣本分別按3:1劃分訓練集和驗證集，各退耕年限得到24個建模樣本和8個驗證樣本，故本研究建模集為96個樣本，驗證集為32個樣本。使用混淆矩陣驗證模型的分類精度，并且使用總精度和Kappa系數對模型精度進行評價，計算公式如下：

式中為混淆矩陣中列的數量；P為混淆矩陣中第行第列的像元數，表示正確分類的個數；P和P分別為第行和第列總樣本個數；代表驗證樣本的總個數。

2 結果與分析

2.1 不同年限退耕土壤的基本理化性質

從表1可知，研究區總樣本土壤有機碳介于2.12～30.11 g/kg之間，變異系數為60.52%，屬于中等變異程度，不同退耕年限土壤有機碳含量依次為RF＞RF＞RF＞RF。含水率最大值為25.39%，最小值為3.77%。飽和導水率介于0.01～10.98 mm/min之間，變異系數為219.61%，屬于強變異程度。黏粒質量百分數介于4.03%～17.81%之間，變異系數為26.46%，屬于中等變異程度。砂粒質量百分數最大值為55.95%，最小值為4.70%，不同退耕年限砂粒質量百分數依次為RF＞RF＞RF＞RF。

表1 土壤樣品基本理化性質Table 1 Basic physical and chemical properties of soil samples

2.2 不同年限退耕土壤的光譜特征

為了分析不同退耕年限土壤光譜曲線的特征，對不同退耕年限土壤光譜曲線求平均值，得到不同退耕年限土壤光譜反射率平均值曲線（圖2a）。從圖2a可以看出，不同退耕年限的土壤光譜曲線形狀基本相似，總體上呈現遞增的趨勢。在可見光波段范圍內，隨著波長的增加，光譜反射率迅速增大，在近紅外波段范圍內，隨著波長的增加，光譜反射率增加的速度逐漸減緩。總體上反射率依次為RF＞RF＞RF＞RF，與SOC含量呈負相關（≤-0.65）。在可見光（400～580 nm）波段則表現為RF＞RF＞RF＞RF，砂粒含量越高，光譜反射率越大。在 610～900 nm 波段范圍內，反射率依次為RF＞RF＞RF＞RF，這可能與土壤中氧化鐵的含量有關，有研究指出，在可見光波段內氧化鐵含量增加會導致反射率明顯下降。本研究對不同退耕年限土壤光譜反射率平均值曲線做去包絡線處理，從而深入分析不同退耕年限土壤的光譜特征（圖2b）。

圖2 不同年限退耕土壤反射光譜曲線及去包絡線Fig.2 Reflectance spectrum curves and continuum removal of soil with different years of conversion

從圖2b中可以看出，在480、900、1 100、1 400、1 900、2 200、2 350 nm處去包絡線都出現了明顯的吸收特征。在400～600 nm波段范圍內，去包絡線值依次為RF＞RF＞RF＞RF，與SOC含量呈正相關，并且表現為吸收谷面積隨著去包絡線值的增大而減小。在800～1 000 nm和1 000～1 200 nm波段范圍內，去包絡線值表現出的規律相同，均為RF＞RF≈RF＞RF，造成RF和RF去包絡線近乎相同的原因可能是RF和RF土壤中氧化鐵含量相近，常認為400～1 100 nm處吸收谷的出現與氧化鐵的存在有著一定的關系。在1 400 nm處的吸收谷是受到土壤表面吸附水、黏土礦物O-H結構晶體等吸收引起的，1 900 nm處的吸收谷主要是由土壤中水分子的O-H官能基在1 900 nm處發生一級倍頻伸縮震動和轉角震動造成的，土壤含水率增加，光譜反射率降低，同時土壤含水率的增加也與黏粒的增加有著一定的關系。對比不同年限退耕土壤飽和導水率可以發現，飽和導水率依次為RF＞RF＞RF＞RF，這與1 900 nm處吸收谷最低點去包絡線值的規律相同。在2 200～2 300 nm附近存在著Al-OH黏土礦物（高嶺石）的吸收帶，從而導致了光譜曲線明顯的吸收谷。已有研究表明，不同粒徑土壤的光譜差異隨波長的增加而增大，因此2 000～2 500 nm處的光譜差異很可能由土壤黏粒、粉粒的差異引起。在2 350 nm處吸收谷的最低點去包絡線值依次為RF＞RF＞RF＞RF，與土壤黏粒含量呈負相關。對土壤理化性質和光譜曲線進行綜合分析，發現不同年限退耕土壤光譜曲線的特征有較大差異，這些差異不僅僅是單一土壤屬性造成的，因此不同年限退耕還林土壤的光譜分類不可以僅依據單一土壤屬性進行。

2.3 不同年限退耕土壤的分類

2.3.1 基于K-means的分類

從聚類結果（表2）可以看出，不同輸入因子構建的聚類模型，聚類精度有很大差距。以CR-PCA作為輸入因子時，模型精度最高，總精度達到81.25%，Kappa系數為0.75。其中RF有1個樣本錯分為RF，RF有1個樣本錯分為RF，RF有3個樣本錯分為RF，僅有1個樣本錯分為RF，并且除RF的聚類精度僅有50.00%外，其余不同年限退耕土壤聚類精度都達到85%以上，說明以CR-PCA為輸入因子，K-means為聚類模型時對不同年限退耕土壤的整體聚類效果較好，對RF的聚類效果則較差。以R-PCA作為輸入因子時，模型精度最差，總精度為46.88%，Kappa系數僅為0.29，僅有RF的聚類精度達到60%以上，其中RF有2個樣本錯分為RF，RF有3個樣本錯分為RF，RF有4個樣本錯分為RF，RF有2個樣本錯分為RF。以RL-PCA、FD-PCA、SCP作為輸入因子時，模型聚類精度較為均衡，從整體上來看，FD-PCA和CR-PCA是K-means模型中相對較好的分類指標，RF、RF以及RF均以CR-PCA的聚類精度為最高，而RF的聚類精度則以FD-PCA為最高。

2.3.2 基于SVM的分類

從表2中可以看出，以R-PCA、RL-PCA為輸入因子時，模型精度較低，總精度和Kappa系數分別為50.00%、46.88%和0.33、0.29，其中RF的分類精度均只有25.00%，為本研究最低的分類精度，并且有3個樣本錯分為RF，2個樣本錯分為RF，1個樣本錯分為RF，說明以SVM作為分類模型，R-PCA和RL-PCA作為輸入因子時，RF的分類效果較差。以FD-PCA作為輸入因子時，模型精度最高，總精度為84.38%，Kappa系數為0.79，不同年限退耕土壤的分類精度都達到了75.00%及以上，RF的分類精度更是達到100%，其中RF有1個樣本錯分為RF，有1個樣本錯分為RF，RF有1個樣本錯分為RF，有1個樣本錯分為RF，RF有1個樣本錯分為RF，說明FD-PCA可以顯著地提高不同年限退耕土壤的分類精度。同K-means分類模型相似，以CR-PCA、SCP為輸入因子時，模型分類精度較為均衡。總體上來看，以SVM作為分類模型時，FD-PCA會顯著提高不同年限退耕土壤的分類精度，CR-PCA、SCP次之。

2.3.3 基于LDA的分類

由表 2可知，分類結果總精度依次為CR-PCA＞FD-PCA=SCP＞R-PCA＞RL-PCA。對比以R-PCA和RL-PCA為輸入因子時不同退耕年限的分類精度可以發現，除RF的分類精度達到80%以上，其余退耕年限土壤的分類精度均小于65%，RF的分類精度只有37.5%，并且RF錯分的樣本大多為RF和RF。以CR-PCA作為輸入因子時，除RF的分類精度較低外，其余退耕年限土壤的分類精度均達到85%以上，具有較好的分類效果，其中RF有1個樣本錯分為RF，RF有1個樣本錯分為RF，1個樣本錯分為RF，RF有1個樣本錯分為RF。以SCP作為輸入因子時，RF和RF分類精度差異較大，其中RF分類精度為100.00%，RF分類精度為62.50%，RF中有1個樣本錯分為RF，2個樣本錯分為RF。以FD-PCA作為輸入因子時，模型分類效果較為均衡。從整體來看，LDA作為分類模型，對RF的分類效果較差；以CR-PCA和FD-PCA為輸入因子時，分類精度有較大提升。

表2 不同模型分類結果Table 2 Classification results of different models

2.4 分類精度比較

以R-PCA、RL-PCA為輸入因子時，精度最高的模型分別是LDA模型和K-means模型，總精度和Kappa系數分別為59.38%、0.46和62.50%、0.50。以FD-PCA、CR-PCA、SCP為輸入因子時，除FD-PCA以SVM模型精度最高，其余均為LDA模型精度最高，總精度分別為84.38%、87.50%、81.25%，Kappa系數分別為0.79、0.83、0.75，一定程度上提高了不同退耕年限土壤光譜分類的精度。在K-means、SVM、LDA模型中，以FD-PCA、CR-PCA為輸入量的模型要略優于以其他三種輸入量的分類模型。從整體上來看，以CR-PCA為輸入因子，基于LDA模型構建的分類模型總精度最高，達到87.50%；以FD-PCA為輸入因子，基于SVM構建的分類模型精度次之，也達到了84.38%；以R-PCA為輸入因子構建的K-means分類模型和以RL-PCA為輸入因子構建的SVM分類模型總精度最差，僅僅只有46.88%；通過對比五種輸入因子在三種模型中的分類總精度可以發現，以CR-PCA和FD-PCA為輸入因子的分類總精度均達到了75%以上，其中CR-PCA總精度最大達到87.50%，FD-PCA總精度最大為84.38%，并且對比不同年限聚類精度，CR-PCA除K-means模型中的RF聚類精度較低外，其余聚類精度均要優于FD-PCA，因此在本研究中CR-PCA是區分不同年限退耕土壤的最優輸入因子。

3 討 論

土壤光譜曲線是土壤屬性的綜合體現。對不同年限退耕土壤的光譜曲線進行分析，發現反射光譜曲線形狀大致相似，反射率依次為RF＞RF＞RF＞RF，與SOC含量呈負相關，與南鋒等研究結果相符，南鋒等對黃土高原煤礦復墾區農田SOC進行預測，發現SOC含量與反射率呈負相關。在400～580 nm波段范圍內，反射率大小與土壤砂粒含量呈正相關，這與張雅梅等的研究結果不符，張雅梅等對土壤質地不同粒級顆粒含量進行統一預測，發現在380～1 100 nm波段范圍內，砂粒含量最高的砂土反射率最低。這可能是SOC含量和砂粒含量共同作用產生的結果。對光譜曲線做去包絡線處理，發現在480和900 nm處，隨著SOC含量的增加，吸收谷變淺，面積減小，這與趙明松等的研究結果相一致，趙明松等對江蘇中部水稻土和潮土的去包絡線進行分析，發現480、900 nm處SOC含量增加導致吸收谷減小。在800～1 200 nm波段范圍內，土壤光譜曲線也表現出一定的差異，這可能由于氧化鐵含量的不同所造成的。去包絡線在1 400、1 900 nm處的吸收谷也表現出顯著差異，造成這一差異的主要原因可能是由于土壤水分的不同。

目前研究多采用高光譜非成像技術、成像技術及高分影像等對土壤類型進行區分，未能考慮諸如退耕還林地不同年限的區分。本文以地面高光譜數據為基礎，實現了退耕還林地不同年限分類研究，可以為高分影像分類時的光譜變換形式和波段選擇提供一定參考，并且后續研究中可以將圖像與光譜結合，探究分類精度更高的模型。本研究以R-PCA、RL-PCA、FD-PCA、CR-PCA和SCP為輸入指標，基于K-means、SVM和LDA構建退耕還林地不同年限分類模型。通過對比不同分類模型的精度，發現以CR-PCA為輸入指標構建的LDA分類模型精度最高。造成這一結果的原因可能是LDA屬于監督學習的線性分類方法，對線性問題有較好的處理效果。比較不同分類模型中輸入指標的精度發現，CR-PCA在K-means、SVM和LDA模型中的分類精度依次為81.25%、78.13%和87.50%，這充分說明對光譜數據進行CR處理可以顯著增強土壤屬性對光譜曲線產生的綜合影響，進而提高模型的精度，其中主要波段范圍為405～595、805～995、1 355～1 495、1 835～2 135、2 145～2 285 nm，后期利用高分影像數據時可以著重考慮選擇這些波段范圍。孟祥添等等采用高分5號影像數據，對原始光譜曲線進行去包絡線、主成分分析等處理，進而區分東北典型黑土區土壤類型，發現CR-PCA分類精度較原始反射率提高了9.15%，說明去包絡線可以增強光譜差異性，與本研究結果一致。

本研究將高光譜技術應用到退耕還林土壤研究中，并實現了退耕還林地不同年限識別，為研究退耕還林土壤屬性變化及影響提供了一種思路和方法。然而本研究構建的模型較為簡單，未構建目前廣泛應用的機器學習模型，會對分類模型的精度有一定的影響，后續將會對分類模型和輸入指標做進一步的研究，為退耕還林對土壤屬性影響的深入研究及退耕還林工作奠定基礎。

4 結 論

為了明確退耕工程對土壤屬性的影響，以及不同年限退耕土壤的光譜特征，實現快速獲取退耕年限的目的，本研究以大寧縣不同年限退耕還林地土壤為研究對象，獲取土壤理化性質及光譜曲線數據，探討不同年限退耕土壤理化性質的變化及對光譜曲線造成的影響，同時構建分類模型，選取最優分類模型及最優輸入因子。主要結論：

1）隨著退耕年限的增長，土壤有機碳含量逐漸增加，砂粒含量先增加后減少，土壤有機碳含量與原始反射率大小呈負相關。

2）不同年限退耕土壤原始光譜曲線形態相似，整體呈遞增的趨勢；進行去包絡線處理后顯著提升光譜曲線的吸收特征，在480 nm處出現了由土壤有機碳引起的吸收谷，900和1 100 nm的吸收谷可能是由于氧化鐵含量不同導致的，1 400和1 900 nm主要是由于土壤水分引起的吸收特征，2 200和2 350 nm處的吸收特征可能是由于不同粒級顆粒含量的不同所導致的。

3）以去包絡線主成分為輸入量，基于線性判別分析構建的分類模型精度最高，達到87.50%，為本研究最優分類模型；以去包絡線主成分為輸入因子的分類模型，總精度均達到了75%以上，最大為87.50%，說明在本研究中去包絡線主成分是區分不同年限退耕土壤的最優輸入因子。