基于無人機高光譜分數階微分的馬鈴薯地上生物量估算

劉 楊 馮海寬 孫 乾 楊福芹 楊貴軍

(1.北京農業信息技術研究中心農業農村部農業遙感機理與定量遙感重點實驗室， 北京 100097；2.山東科技大學測繪科學與工程學院， 青島 266590； 3.國家農業信息化工程技術研究中心， 北京 100097；4.河南工程學院土木工程學院， 鄭州 451191； 5.北京市農業物聯網工程技術研究中心， 北京 100097)

0 引言

及時掌握馬鈴薯的生長狀況對于指導田間管理和優化種植格局、挖掘生產潛力具有重要意義。地上生物量(AGB)是反映作物生長狀況的重要指標，AGB的變化能夠直接表征作物有機物的積累能力，反映作物生長狀況和營養狀態[1-2]。快速準確地獲取AGB信息對于監測馬鈴薯光合作用能力和生長狀態至關重要[3-6]。傳統的人工測量手段獲取AGB信息準確，但破壞性極強，且難以滿足大規模種植的監測需求。高光譜遙感技術具有較高的波譜分辨率和較強的波段連續性，是對地觀測的重大突破[7-10]。由于作物冠層對太陽光的吸收和反射形成特有的光譜曲線，故可通過分析作物的高光譜反射特性實現作物AGB的快速無損估測[11-12]。

無人機高光譜遙感具有機動性強、成本低、可云下獲取影像等優勢，已成為當前農業領域所關注的焦點。國內外學者基于無人機高光譜信息在作物地上生物量估算研究方面取得了一定的成果[13-19]。

現有研究多基于原始冠層光譜特征或利用其構建的植被指數與地上生物量建立關系，構建作物地上生物量估算模型。光譜微分技術如分數階微分能夠細化光譜信息，消除環境背景的影響，深度挖掘光譜中潛在的信息。近年來，利用分數階微分技術在各個領域得到了廣泛研究與應用[20-23]。然而，分數階微分技術在作物營養監測方面的研究還較少。

本文對波段454～950 nm的馬鈴薯冠層無人機高光譜數據進行0～2階(間隔0.2)微分處理，得到每個生育期11種冠層分數階微分光譜數據。通過相關性分析，挑選出各生育期與馬鈴薯AGB相關性較好的各階分數階微分光譜，計算其相關系數絕對值，優選不同階微分波段，構建馬鈴薯AGB的多元線性回歸模型、隨機森林模型和人工神經網絡模型，并進行模型精度驗證，最后篩選出最優估算模型。

1 試驗設計與數據獲取

1.1 試驗設計

于2019年3—7月在北京市昌平區小湯山鎮國家精準農業研究示范基地的馬鈴薯試驗田開展試驗，位于北緯40°10′34″，東經116°26′39″，平均海拔為36 m，氣候類型為暖溫帶半濕潤大陸性季風氣候。如圖1所示，試驗區共設密度試驗(N區)、氮素試驗(S區)、鉀肥試驗(K區)3個試驗區，采用完全隨機試驗設計，N區和S區均設2個試驗品種，分別為中薯5(P1)和中薯3(P2)，均為早熟品種，K區設1個試驗品種，為中薯3。密度試驗設3個水平：60 000株/hm2(T1)、72 000株/hm2(T2)、84 000株/hm2(T3)，6個處理，每個處理重復3次，共18個試驗小區；氮素試驗設4個水平(以尿素計量)：0 kg/hm2(N0)、244.65 kg/hm2(N1)、489.15 kg/hm2(N2,正常處理，15 kg純氮)、733.5 kg/hm2(N3)，8個處理，每個處理重復3次，共24個試驗小區；鉀肥試驗設3個水平：0 kg/hm2(K0)、970.5 kg/hm2(K1，N區和S區均為K1處理)、1 941 kg/hm2(K2)，2個處理，每個處理重復3次，共6個試驗小區；試驗小區總計48個，試驗小區面積為6.5 m×5 m。

圖1 馬鈴薯試驗位置和試驗設計Fig.1 Experimental area and design of potato

1.2 地面數據獲取

分別獲取馬鈴薯現蕾期(2019年5月13日)、塊莖形成期(2019年5月28日)、塊莖增長期(2019年6月10日)、淀粉積累期(2019年6月20日)、成熟期(2019年7月3日)5個關鍵時期的地上生物量數據。馬鈴薯地上生物量通過收獲法獲取，在每個小區中選取代表性的3顆植株，將其莖葉分離，隨后用清水洗凈，105℃殺青，80℃干燥48 h以上，直到質量恒定再進行稱量。將植株莖、葉的干質量求和得到樣本干質量，最后通過樣本干質量與種植密度得到每個小區的馬鈴薯地上生物量[24-25]。

1.3 無人機高光譜數據獲取及處理

試驗利用無人機搭載的德國Cubert公司生產的UHD185型機載成像光譜儀，該光譜儀尺寸為195 nm×67 nm×60 nm，質量470 g，波譜范圍為450～950 nm，共有125個光譜通道，采樣間隔4 nm，光譜分辨率8 nm，數字分辨率12位。選擇晴朗無云的天氣進行拍攝，此時太陽光照強度穩定。飛行時間10:30—14:00，飛行高度為50 m。

無人機獲取的高光譜數據處理主要包括輻射校正、影像拼接、影像融合和光譜提取4部分。根據UHD185型成像光譜儀中心波長和波長半幅寬在Matlab環境設計的輻射定標系統，完成由影像DN值到地表反射率的輻射校正。利用Agisoft PhotoScan 軟件進行高光譜影像拼接，變為cue格式數據，再提取子波段為jpg格式，最后再對各個拼接的子波段影像進行合并生成馬鈴薯試驗田的高光譜影像數據。基于Arcmap10.2軟件，繪制馬鈴薯每個小區矢量數據，基于IDL語言統計每個小區平均光譜作為馬鈴薯冠層光譜，得到高光譜反射率數據。

1.4 數據分析方法

分數階微分是從整數階微分發展而來的，本質上是任意階光譜反射率的斜率。目前，比較流行的分數階微分定義為Grunwald-Letnikov、Riemann-Liouville和Caputo 3種形式[15]。3種定義中，Grunwald-Letnikov因形式簡單而被廣泛運用，故本文運用此形式完成對馬鈴薯冠層光譜信息的深度挖掘，該微分形式為

(1)

式中Γ——Gamma函數λ——對應波長

n——微分上下限之差

a——任意階數f(λ)——光譜

當a=0、1或2時，f(λ)為原始光譜、一階微分光譜或二階微分光譜；當a是小數時，則為Grunwald-Letnikov微分分數階形式。

1.5 精度評價

本文對每個生育期分別挑選2/3樣本數據作為建模集，1/3樣本數據作為驗證集，以此構建馬鈴薯AGB估算模型。采用決定系數(Coefficient of determination,R2)、均 方 根 誤 差(Root mean square error, RMSE)、標準均方根誤差(Normalized root mean square error, NRMSE)評價模型的精度。R2越接近于1，RMSE和NRMSE越低，其估測模型的精度就越高。

2 結果與分析

2.1 原始冠層光譜與生物量相關性

由圖2可知，馬鈴薯冠層光譜曲線呈現典型的綠色植物反射特征，在可見光波段范圍內，存在2個吸收帶，分別位于480 nm和670 nm附近，在550 nm附近形成小的反射峰，這一特征變化主要由于馬鈴薯葉片中葉綠素對綠光反射作用強，對藍光和紅光吸收強造成的。在波段674～780 nm之間馬鈴薯冠層反射率急增，形成典型綠色植被特有的“紅邊”特征。

圖2 不同試驗的馬鈴薯冠層光譜曲線Fig.2 Potato canopy spectral curves of different experiments

為篩選出與馬鈴薯AGB相關性較好的原始冠

層光譜波段，通過相關性分析，得到各個生育期原始冠層光譜與AGB的相關性，結果如圖3所示。從圖3可以看出：①現蕾期，冠層原始光譜在波段454～718 nm范圍內與馬鈴薯AGB呈極顯著負相關(P<0.01)，在波段746～914 nm范圍內與AGB呈極顯著正相關。由于與AGB相關的光譜波段主要是可見光波段，故選取相關性較大的波長478、674 nm，其相關系數分別為-0.64和-0.66。②塊莖形成期，冠層原始光譜分別在波段558～698 nm、726～950 nm與AGB呈極顯著負、正相關。相關性較高的波長為766、778 nm，其相關系數分別為0.63和0.63。③塊莖增長期，冠層原始光譜分別在波段454～702 nm、718～950 nm與AGB呈極顯著負、正相關。相關性較大的波長為610、770 nm，其相關系數分別為-0.72和0.71。④淀粉積累期，冠層原始光譜分別在波段454～710 nm、714～950 nm與AGB呈極顯著負、正相關。相關性較高的波長為482、618、746 nm，其相關系數分別為-0.70、-0.67和0.76。⑤成熟期，冠層原始光譜在波段914～950 nm與AGB呈極顯著負相關，未見波段與AGB呈極顯著正相關。在可見光波段范圍內選擇相關系數絕對值較大的波長550、694 nm，其相關系數為-0.33和0.25。

圖3 馬鈴薯不同生育期原始冠層光譜與AGB相關系數Fig.3 Correlation between original canopy spectrum of potato and AGB at different growth stages

2.2 分數階微分光譜與生物量相關性分析

為篩選出每個生育期與馬鈴薯AGB相關性較好的分數階微分光譜波段，通過相關性分析，得到不同生育期馬鈴薯冠層分數階微分光譜與AGB的相關性，0.4階、0.8階、1.2階、1.6階和2.0階的相關系數變化如圖4所示。

根據馬鈴薯不同生育期冠層分數階微分光譜與AGB的相關性，得到不同分數階下微分光譜與AGB相關系數絕對值最大值，繪制不同生育期下不同階微分光譜與AGB相關系數絕對值最大值的折線圖，結果如圖5所示。分析圖5可知，不同生育期馬鈴薯AGB與分數階微分光譜間的相關系數絕對值最大值出現的階數不同。現蕾期，相關系數絕對值最大值在0.8階微分(470 nm)；塊莖形成期，相關系數絕對值最大值在1.8階微分(710 nm)，塊莖增長期和淀粉積累期，相關系數絕對值最大值都在1.6階微分(718、722、766 nm)；成熟期，相關系數絕對值最大值在1.0階微分(622 nm)。因此表明常用的“綠邊(502～554 nm)”和“紅邊(678～758 nm)”雖對作物營養狀況能夠很好地進行監測，但僅采用二者范圍內的波段進行估算AGB時，并不能夠深度挖掘光譜中潛在的信息，導致AGB估算時產生飽和現象，然而利用分數階微分光譜能夠凸顯波段中隱藏的信息，提高光譜信息的敏感度，增強“綠邊”和“紅邊”位置對AGB估算的精度。

圖4 馬鈴薯不同生育期的分數階微分光譜與AGB的相關系數Fig.4 Correlation between fractional differential spectrum of potato at different periods and AGB

圖5 不同生育期相關系數絕對值最大值出現的階數Fig.5 Order of maximum absolute value of correlation coefficient in different growth periods

2.3 敏感波段確定

為了盡可能避免各個分數階微分光譜之間的共線性，將與AGB相關系數絕對值按從大到小排列，挑選出各個生育期與馬鈴薯AGB相關的前9個分數階微分冠層光譜波段，并將各生育期的相關系數絕對值繪制成矩陣圖，結果如圖6所示。由圖可知，現蕾期，0.2階微分在波長666 nm處，0.4階微分在波長634 nm處，0.6階微分在波長470、578、594 nm處，0.8階微分在波長458、470 nm處，1.2階微分在波長790 nm處，1.6階微分在波長706 nm處，與馬鈴薯AGB達到0.01顯著水平，相關系數絕對值在0.67～0.72之間；塊莖形成期，1.2階微分在波長714 nm和718 nm處、1.4階微分在波長710 nm和718 nm處、1.6階和1.8階微分都在波長710 nm和714 nm處、2.0階微分在波長710 nm處，與馬鈴薯AGB達到0.01顯著水平，相關系數絕對值在0.72～0.76之間；塊莖增長期，0.8階微分在波長766 nm處、1.0階微分在波長742 nm處、1.2階微分在波長730 nm和894 nm處、1.4階微分在波長726 nm處、1.6和1.8階微分都在波長718 nm和722 nm處，與馬鈴薯AGB達到0.01顯著水平，相關系數絕對值均在0.76以上；淀粉積累期，0階微分在波長746 nm處、0.2階微分在波長730 nm和742 nm處、0.4階微分在波長726 nm處、0.6階微分在波長730 nm處、0.8階微分在波長726 nm處、1.0階和1.2階微分都在波長722 nm處、1.6階微分在766 nm波長處，與馬鈴薯AGB達到0.01顯著水平，相關系數絕對值在0.75～0.78之間；成熟期，0.2階微分和0.4階微分都在922 nm波長處，0.6階微分在波長686 nm和690 nm處，0.8階微分在波長622、682、686 nm處，1.0階微分和1.2階微分都在波長622 nm處，與馬鈴薯AGB達到0.01顯著水平，相關系數絕對值在0.50～0.59之間。綜合每個生育期挑選的前9個微分光譜波長可知，僅采用整數階微分并不能充分表達與AGB之間的聯系，使得冠層光譜中隱藏的信息被遺漏，而通過分數階微分能夠細化光譜信息，深度挖掘出與AGB相關的有效信息，使得分數階微分光譜較整數階微分光譜與AGB的相關性更高，提高AGB估算模型精度。

圖6 不同生育期馬鈴薯AGB與分數階微分光譜相關系數矩陣Fig.6 Correlation coefficient matrix of potato AGB and fractional differential spectrum in different periods

2.4 基于分數階微分光譜的馬鈴薯生物量最優模型篩選

以馬鈴薯地上生物量(AGB)為因變量，由0～2階分數階微分(0.2間隔)確定的前9個微分波段為自變量，建立分數階微分冠層光譜與馬鈴薯AGB的多元線性回歸(MLR)模型、隨機森林(RF)模型和人工神經網絡(ANN)模型，并驗證模型精度，進而挑選出各個生育期最優估算模型，結果見表1。由表1可知，利用MLR、RF和ANN方法估算馬鈴薯AGB具有較高的精度，從現蕾期到塊莖增長期這3個生育期估算效果逐漸變優，驗證R2也是不斷增大，RMSE和NRMSE逐漸降低，從淀粉積累期到成熟期建模和驗證R2逐漸降低，RMSE和NRMSE逐漸增加，每種方法均在塊莖形成期估算效果最佳。每個生育期的AGB估算，通過MLR構建的模型效果最優，其次為RF模型，而ANN模型效果最差，且采用3種方法估算AGB得到的模型均是驗證效果要優于相應的建模效果。塊莖增長期，利用3種方法建模R2分別為0.76、0.72和0.61，RMSE分別為143.36、135.06、160.44 kg/hm2，NRMSE分別為16.67%、17.63%、19.58%；驗證R2分別為0.84、0.74和0.69，驗證RMSE分別為94.72、101.44、159.09 kg/hm2，驗證NRMSE分別為11.61%、14.02%和17.76%。

為了弄清不同因素對模型的作用，以不同品種、密度和施肥狀況下的馬鈴薯AGB為因變量，采用上述3種方法構建AGB估算模型，并計算各個模型的精度評價指標，其結果見表2～4。綜合表2～4可知，利用MLR、RF和ANN方法估算3種狀況下的馬鈴薯AGB，也是從現蕾期到塊莖增長期這3個生長期估算效果逐漸變優，驗證R2也是不斷增大，RMSE和NRMSE逐漸降低，從淀粉積累期到成熟期建模和驗證R2逐漸降低，RMSE和NRMSE逐漸增加，每種方法均在塊莖形成期估算效果最佳。每個生育期不同狀況的AGB估算，也是通過MLR構建的模型效果最好，其次為RF模型，而ANN模型效果最差。每種狀況下AGB估算，采用3種方法得到的模型均是驗證效果要優于相應的建模效果，其中以不同品種的AGB模型，每種方法以中薯5(P1)構建模型的建模和驗證精度要優于相應的以中薯3(P2)所構建的模型精度。綜上可知，整體上得到的模型效果與不同狀況下的模型效果基本一致，說明品種、密度和施肥對AGB估算具有等效的影響力。

表1 不同生育期的馬鈴薯AGB估測精度對比Tab.1 Comparison of accuracy of potato AGB estimation at different growth stages

表2 不同品種的馬鈴薯AGB估算精度對比Tab.2 Comparison of accuracy of AGB estimation for different varieties of potatoes

表3 不同密度的馬鈴薯AGB估算精度對比Tab.3 Comparison of accuracy of AGB estimation for different densities of potatoes

表4 不同施肥的馬鈴薯AGB估算精度對比Tab.4 Comparison of accuracy of potato AGB estimation under different fertilizations

綜合建模和驗證集可知，每個生育期建模R2低于驗證R2，所得R2不僅與精度有關，也與樣本數量有關，整體、不同品種、不同密度和不同施肥條件下建模樣本分別為32個、(14個、18個)、12個和20個，驗證樣本分別為16個、(7個、9個)、6個和10個，驗證樣本的數量小于建模數量，所以驗證R2高于建模的R2；驗證模型的RMSE和建模的RMSE接近，且二者得到的NRMSE基本處于20%以內，說明模型的穩定性較好，預測能力較高。

每個生育期采用3種方法都在塊莖增長期表現效果最好，主要因為馬鈴薯從現蕾期開始生長，最初為營養生長和生殖生長，體現在莖節迅速伸長，葉片快速擴大。到了塊莖增長期，莖節和葉片等營養器官生長最為旺盛，此時植被覆蓋度最大，冠層光譜信息的提取不易受到地面土壤的影響，使得冠層光譜信息能夠充分表達和AGB之間的聯系，進而構建的模型效果最優。而淀粉積累期后，地上部有機物不斷向地下輸送，外加連續多天大雨，造成地上部葉片迅速枯黃脫落，植被覆蓋度明顯降低，冠層光譜信息的提取容易受到地面土壤的影響，使得冠層光譜信息與AGB的聯系變差，進而構建的模型效果不好。

每個生育期利用MLR方法構建馬鈴薯AGB估算模型的建模精度和驗證模型的穩定性優于RF和ANN模型，主要因為機器學習適用于處理較大數據集，而對較小數據集表現優勢不明顯，本文建模樣本和驗證樣本都屬于較小數據集，所以RF和ANN的表現能力較差一些。而ANN模型估算精度和驗證模型的穩定性最差，可能是利用ANN方法構建模型時，在訓練過程中不斷地反復學習，過多地學習了細節，造成估算能力不佳。每個生育期利用3種方法估算馬鈴薯AGB時，每種方法的建模精度由高到低依次為塊莖增長期、塊莖形成期、淀粉積累期、現蕾期、成熟期。

3 結論

(1)以馬鈴薯為研究對象，利用無人機冠層高光譜數據，采用分數階微分方法以0～2階(間隔0.2)進行光譜微分處理，篩選出前9個分數階微分光譜，結合地上實測的48個生物量數據，建立了估算馬鈴薯生物量的MLR、RF和ANN模型。

(2)在不同生育期，與AGB相關性較高的原始冠層光譜波段不同。現蕾期，在可見光范圍內，與AGB極顯著相關的波長為478、674 nm，其相關系數分別為-0.64、-0.66；塊莖形成期，與AGB極顯著相關的波長為766、778 nm，其相關系數均為0.63；塊莖增長期，與AGB極顯著相關的波長為610、770 nm，其相關系數分別為-0.72、0.71；淀粉積累期，與AGB極顯著相關的波長為482、618、746 nm，其相關系數分別為-0.70、-0.67和0.76；成熟期，冠層原始光譜在可見光范圍內未存在極顯著正負相關的波長，而相關性較好的波長為550、694 nm，其相關系數為-0.33、0.25。

(3)不同生育期分數階微分光譜與AGB相關系數絕對值最大值出現的階數不同。現蕾期，相關系數絕對值最大值在0.8階微分(470 nm)，為0.72；塊莖形成期，相關系數絕對值最大值在1.8階微分(710 nm)，為0.76；塊莖增長期和淀粉積累期，相關系數絕對值最大值都在1.6階微分(718、722、766 nm)，分別為0.77和0.78；成熟期，相關系數絕對值最大值在1.0階微分(622 nm)，為0.59。通過比較馬鈴薯成熟期的原始光譜和微分光譜與AGB的相關性可知，分數階微分光譜更能挖掘出光譜數據中隱藏的信息。

(4)由每個生育期挑選的前9個微分光譜波段可知，相較于整數階微分，分數階微分與AGB的相關性更高。基于分數階微分光譜的AGB估算模型，不同生育期利用MLR方法，以整體、不同品種、不同密度和不同施肥條件下的AGB為自變量，所得模型效果最好，其次是RF模型，而ANN模型效果最差。不同生育期利用3種方法估算馬鈴薯AGB時，其建模精度由高到低依次為塊莖增長期、塊莖形成期、淀粉積累期、現蕾期、成熟期。