基于無人機高光譜影像的冬小麥全蝕病監(jiān)測模型研究

郭 偉 朱耀輝 王慧芳 張 娟 董 萍 喬紅波

(1.河南農業(yè)大學信息與管理科學學院， 鄭州 450002； 2.北京市氣象局北京氣候中心， 北京 100089)

0 引言

小麥全蝕病是世界各小麥產區(qū)危害十分嚴重的一種毀滅性病害，全蝕病的發(fā)生會引起植株成簇或大片枯死，嚴重時減產一半以上，甚至絕收，造成嚴重的產量損失，因此防治小麥全蝕病對保障國家糧食安全至關重要[1-4]。全蝕病為土傳病害，一旦發(fā)病，次年發(fā)病概率較高，且預防來年小麥全蝕病的拌種劑價格昂貴。針對這一特點，在空間上及時、精確地監(jiān)測當年小麥全蝕病的發(fā)病區(qū)域，估測病害嚴重度，可預測下一年的易發(fā)病區(qū)域，精準施藥，降低生產成本，從而指導人們防治全蝕病病害。

小麥全蝕病由根部感染全蝕病菌引起，幼苗期植株矮小，灌漿期 “白穗”是其典型病狀。采用遙感技術可以快速、準確獲取其染病區(qū)域的空間信息。目前，基于傳統遙感技術的作物病蟲害監(jiān)測已經廣泛應用于作物病害監(jiān)測。但是在獲取農作物冠層光譜時，傳統的非成像近地光譜儀每次只能獲取單點光譜數據，而近地成像光譜儀受到平臺高度的限制，獲取的影像范圍非常有限，且在大田環(huán)境中難以獲取正射影像；遙感衛(wèi)星雖可以快速獲取大面積農田的影像，但受衛(wèi)星過境時間、天氣等條件限制較多，且成本較高[5]；近年來，無人機遙感技術以其數據測量時間靈活機動、時空分辨率高、觀測范圍大等優(yōu)勢，在農田生態(tài)環(huán)境信息監(jiān)測與獲取中開始應用[6-9]，但目前農業(yè)無人機遙感多以普通數碼相機和多光譜相機為主要傳感器，所獲取的影像數據波段少、光譜信息有限，這些不足限制了遙感技術在農作物病蟲害監(jiān)測中的應用。低空無人機搭載高分辨率成像光譜儀對小麥病蟲害進行空間監(jiān)測可以很好地彌補以上不足。近年來，國內外已經陸續(xù)出現了利用無人機高光譜成像儀進行農田監(jiān)測的研究[10-12]。初步評價了無人機成像高光譜數據在農業(yè)應用方面的潛力，并利用無人機高光譜技術對作物長勢及生理生化參數進行估算。但在作物病蟲害空間監(jiān)測方面的研究鮮有報道。

本研究以低空無人機作為遙感平臺，使用高光譜分辨率的新型成像光譜儀，以患全蝕病冬小麥為研究對象，結合地面調查數據和近地高光譜數據，客觀評價無人機高光譜遙感數據質量，運用高光譜技術，綜合分析不同光譜指數和最小二乘法在冬小麥全蝕病病情指數估測的能力，并通過低空無人機獲取的高光譜成像數據實現冬小麥全蝕病病情指數的空間反演分布，以期為基于無人機高光譜遙感的作物病蟲害空間監(jiān)測積累經驗。

1 材料與方法

1.1 試驗地概況

試驗在河南省中南部的漯河市郾城區(qū)裴城鎮(zhèn)(33°43.1′N，113°49′E)進行，試驗地長112 m、寬26 m，如圖1所示。試驗冬小麥品種為矮抗58，于2016年10月12日播種，冬小麥生長期間的栽培措施同一般高產田。

1.2 數據獲取與處理

1.2.1病情指數調查

圖1 研究區(qū)示意圖Fig.1 Schematic of study area

數據獲取時間為2017年5月21日，此時為冬小麥灌漿期，患全蝕病冬小麥的“白穗”特征最為明顯。依據病害程度的不同，隨機選取40個樣點，于每個樣點西南方向45°處插入標記板并使用差分GPS記錄位置信息。獲取樣點冠層光譜數據和無人機高光譜成像數據后，以每一個光譜測定樣點右前方1 m×1 m處作為病情指數調查樣方，按照GB/T 17980.109—2004規(guī)定，每個樣方按照傳統的五點取樣方法選取小麥樣本，每個樣方取1 m長雙行冬小麥調查總穗數及白穗數。取回樣品由植保專業(yè)實驗員在實驗室環(huán)境下提取病害組織置于培養(yǎng)基培養(yǎng)，純化病原菌采用柯赫氏法則確定病原菌，鑒定冬小麥全蝕病病株。清水沖洗樣品根系，在白色背景下調查根系發(fā)病情況，根據每個植株根部受侵染面積的百分比計算病情級數，病情級數分類標準如表1所示。

表1 冬小麥全蝕病病情級數分類標準Tab.1 Classification criteria for winter wheat take-all grade

白穗率從冠層角度反映了冬小麥發(fā)病區(qū)域和面積，病情級數則是從根系角度反映了冬小麥患全蝕病的嚴重程度。病情指數是全面考慮發(fā)病率與嚴重度的綜合指標，作為本文評價冬小麥全蝕病患病程度的綜合指標。其計算公式為

(1)

(2)

式中WHR——白穗率DWH——枯白穗數

TWE——總穗數AGWH——各級白穗數

TWH——總白穗數DG——病情級數

HDG——最高病情級數

DI——病情指數

1.2.2無人機高光譜成像數據

無人機遙感平臺如圖2所示，包含八旋翼無人機、飛行控制系統、成像高光譜儀、微型單板計算機Pokini Z、無線遙感系統、地面控制系統和數據處理系統。其中八旋翼無人機為天途公司生產的AZUP-T8，單臂長60 cm，機身凈質量7 kg，可載質量10 kg，續(xù)航時間為30 min；成像光譜儀為Cubert UHD185，拍攝方式為全畫幅、非掃描、實時成像，光譜范圍為450～950 nm，波段數125，光譜采樣間隔4 nm，光譜分辨率8 nm。無人機遙感作業(yè)在天氣晴朗無云、風速較小時進行，航速6 m/s，航高50 m，航向重疊度80%，旁向重疊度60%，地面控制站通過無線網絡遠程遙控Pokini Z運行，并將拍攝數據存儲在Pokini Z中。

UHD185型機載高速成像光譜儀獲取的數據包含：空間分辨率21 cm的高光譜立方體影像；空間分辨率1 cm的全色jpg圖像。UHD185遙感數據預處理包括兩部分：高光譜立方體影像的拼接和試驗區(qū)域樣點冬小麥冠層平均光譜的提取。采用德國Cubert公司研發(fā)的Cube-Pilot軟件與Agisoft LLC公司研發(fā)的Agisoft PhotoScan軟件進行高光譜影像拼接，高光譜影像拼接流程如圖3所示。

圖2 無人機高光譜設備與遙感影像采集Fig.2 UAV hyperspectral equipment and remote sensing image acquisition

圖3 高光譜影像的拼接流程圖Fig.3 Splicing process of hyperspectral image

1.2.3冠層光譜數據

地面數據采集工作與無人機空中作業(yè)同步開展。在進行冬小麥全蝕病病情指數調查前，先對樣點進行光譜測量。光譜測量采用美國ASD HandHeld手持地物光譜儀，光譜范圍為 325～1 075 nm，光譜采樣間隔為1 nm，光譜分辨率為3 nm，每次采集目標光譜前均進行參考板校正，傳感器探頭垂直向下，且探頭高度始終距離地面1.3 m(此時小麥平均株高為0.78 m)，探頭視場角25°，地面視場范圍直徑為50 cm，每個樣點均隨機采集10次光譜數據，取其平均值作為該采樣點的光譜測量值，并用標記板標記樣點序號，同時用差分GPS記錄樣點經緯度坐標。

以ASD Viewspec-pro軟件處理光譜數據。為了更好地反映不同病害程度在冠層光譜上的響應特點，依據前人研究[13-15]根據不同DI值將其劃分為4個病害等級，分別是正常(0)、輕微(0～30%)、中等(30%～50%)和嚴重(>50%)。

1.3 方法

1.3.1光譜指數的構建

為了更好地反映病情指數與光譜信息之間的關系，本研究從拼接后的UHD185高光譜影像上提取每個樣點的平均光譜數據，在450～950 nm范圍內將任意兩個光譜波段構建光譜指數，以便更好地利用光譜信息。光譜指數是不同光譜波段的植被-土壤系統的反射率因子以一定的代數形式組合而成的一種參數[16],依據小麥患全蝕病后的外部及生理特征和前人研究進展[17-18]選取能夠較好反映植物結構或者葉面積指數、生物量等生理活性參數的歸一化差異光譜指數(NDSI)、差值光譜指數(DSI)和比值光譜指數(RSI)，定義為

(3)

DSI(Ri,Rj)=Ri-Rj

(4)

(5)

式中Ri——450～950 nm波段內的第i個波長對應的光譜反射率

Rj——450～950 nm波段內的第j個波長對應的光譜反射率

實現工具為Matlab 2015b。

通過任意兩波段構建DVI、NDSI和RSI與病情指數的決定系數等勢圖，尋找具有更高預測能力的光譜指數和敏感波段。

1.3.2偏最小二乘回歸

偏最小二乘回歸(Partial least squares regression, PLSR)方法已經被廣泛應用到植被生理參數估測中。其在多元線性回歸方法的基礎上，結合了主成分分析(PCA)和典型相關分析(CCA)的思想，這樣自變量之間相互獨立且能夠最大程度反映與因變量之間的差異，解決了多元線性回歸分析中常常遇到的小樣本以及多重共線性問題。本研究將每個樣點的病情指數(DI)作為因變量，構建的3種光譜指數作為自變量，建立PLSR預測模型，預測病情指數。

1.3.3統計分析

1.3.3.1決定系數

決定系數R2為評價UHD185和重采樣后ASD的光譜反射率的相關性，以及各個樣點DI值與光譜指數構建模型的相關性評價標準。R2越大，表示參與評價的變量間相關性越好，參考價值越大；R2越小，表示兩者相關性較差，參考價值越低。

通過Matlab軟件編輯代碼計算462～874 nm(4～102波段)范圍內任意兩波長構建的差值光譜指數(DSI)、歸一化差異光譜指數(NDSI)和比值光譜指數(RSI)，并構建小麥全蝕病病情指數的R2等勢圖，選出線性關系最優(yōu)的波段組合，分析對比具有較高預測能力的光譜指數。

1.3.3.2均方根誤差和平均絕對誤差

均方根誤差(Root mean square error, RMSE)和平均絕對誤差(Mean absolute error, MAE)被用來作為評價估算病情指數值與實測病情指數值的標準。

對于病情指數，RMSE對較大誤差非常敏感，RMSE越小，表示對病情指數的估算精度越高。MAE將各次測量值的絕對誤差取絕對值后再求平均值，平均絕對誤差由于離差被絕對值化了，平均誤差可以更好地反映預測值誤差的實際情況，同時使用RMSE和MAE有助于更好地分析誤差，避免個別誤差造成RMSE過大。

2 結果與分析

2.1 UHD185影像可靠性評價

UHD185成像高光譜數據由高光譜影像與全色影像融合得到，在利用其進行冬小麥全蝕病監(jiān)測前，評價其光譜數據的可靠性是必要的。與無人機空中拍攝成像高光譜影像同步作業(yè)的ASD地物光譜儀在農業(yè)遙感中已被廣泛使用，其光譜信息往往作為作物長勢和病蟲害監(jiān)測的重要依據[17-19]，因此本文以ASD光譜數據為標準，在每個樣點西南方向插入標記板并使用差分GPS記錄位置信息。利用地面布設標記板來確保采集的ASD與UHD185光譜數據為相同樣點測量區(qū)域，通過對比ASD與UHD185的光譜反射率，對UHD185數據進行可靠性評價。

圖5 不同病害等級的UHD185與重采樣ASD光譜反射率的相關性Fig.5 Correlation between UHD185 and resampled ASD spectral reflectance at different disease grades

首先，將ASD光譜數據重采樣為UHD185波段，分別計算不同病害等級的UHD185與ASD光譜數據的平均光譜反射率，如圖4所示。UHD185光譜曲線與ASD光譜曲線的變化趨勢總體一致。UHD185在550 nm附近出現“綠峰”特征，在680 nm附近出現“紅谷”特征，與ASD光譜曲線相符。原因在于冬小麥冠層光譜反射率隨葉片中葉綠素含量的增加而減少，由于葉綠素水平的升高，相應葉片光合作用增強，對紅藍光吸收能力增強，相應的光譜反射率減少[20]；在680～750 nm的紅邊區(qū)域，兩者的光譜反射率急劇上升，此波段內光譜曲線對應斜率最大值的波長是紅邊位置REP；在770～910 nm的近紅外區(qū)域，相同冬小麥全蝕病病害等級的冠層光譜曲線近似一致，呈現近紅外高反射平臺，不同病害等級間ASD和UHD185光譜反射率由大到小為正常、輕微、中等、嚴重，原因是冬小麥白穗程度越高，冬小麥冠層葉片中水分越少，葉綠素含量、氮含量和葉面積指數減少，導致冠層光譜在近紅外波段的反射率也隨之降低，BONHAM-CARTER[21]的研究驗證了近紅外波段反射率與氮素水平正相關，與本文結果相似；在910 nm之后，ASD和UHD185光譜反射率逐漸降低，且UHD185光譜曲線波動性較大。

圖4 不同病害等級UHD185光譜曲線與重采樣的ASD光譜曲線對比Fig.4 Comparison between UHD185 spectral curves and resampled ASD spectral curves at different disease grades

然后，計算冬小麥冠層的重采樣ASD和UHD185光譜反射率在462～874 nm(4～102波段)波段內的線性關系，結果表明，兩者決定系數R2均在0.97以上，如圖5所示。綜上兩方面的分析得出，在462～874 nm(4～102波段)范圍的UHD185光譜數據最為可靠，可利用其對不同全蝕病冬小麥病害等級進行精準分類。

2.2 基于成像高光譜數據的冬小麥全蝕病病情指數預測模型

2.2.1估測冬小麥全蝕病病情指數的最優(yōu)光譜指數

圖6為任意兩波長組合構成的DSI(Ri,Rj)、NDSI(Ri,Rj)和RSI(Ri,Rj)與病情指數的決定系數等勢圖。

圖6 任意兩波長組合構成的DSI(Ri,Rj)、NDSI(Ri,Rj)和RSI(Ri,Rj) 與病情指數的決定系數等勢圖Fig.6 Contours of R2 between disease index and DSI(Ri,Rj), NDSI(Ri,Rj) and RSI(Ri,Rj)

從圖6的R2等勢圖中可看出估測冬小麥全蝕病病情指數的光譜指數最優(yōu)波段組合。對于DSI(Ri,Rj)(圖6a)，R2>0.5的波段區(qū)域為730～874 nm與462～760 nm的波段組合，其中DSI(R818,R534)波長組合的相關性最好，決定系數R2達到了0.860 5；對于NDSI(圖6b)，R2>0.5的波段區(qū)域為510～710 nm與740～850 nm、550～660 nm與550～580 nm，其中NDSI(R530,R758)波長組合的相關性最好，決定系數R2達到了0.653 6；對于RSI(圖6c)，R2>0.5的波段區(qū)域為510～710 nm與740～850 nm、550～660 nm與550～580 nm，其中RSI(R758,R530)波長組合的相關性最好，決定系數R2達到了0.659 0，如表2所示。

表2 特征波長組合的光譜指數與病情指數相關性Tab.2 Correlation between spectral index and disease index of characteristic band combination

從波段位置看，DSI、NDSI和RSI所優(yōu)選出的波段距550 nm“綠峰”非常近，說明其構建的DSI(R818,R534)、RSI(R530,R758)和NDSI(R758,R530)與病情指數的相關性受葉片色素含量影響較為明顯，這是由患全蝕病冬小麥的白穗特征導致；整體而言，DSI光譜指數與冬小麥全蝕病病情指數相關性最好，DSI與全蝕病冬小麥病害指數的敏感性波段相較于NDSI、RSI更寬，冬小麥全蝕病病情指數與矩陣聯立形式任意波段組合的相關性也最高，驗證了DSI(R818,R534)在精準監(jiān)測冬小麥全蝕病病害的能力上優(yōu)于NDSI(R758,R530)、RSI(R530,R758)。

2.2.2基于光譜指數的冬小麥全蝕病DI模型構建及檢驗

基于本研究中UHD185高光譜影像提取的40個樣點數據，其中20個樣點作為建模集，另外20個樣點作為驗證集，以冬小麥全蝕病DI為因變量，分別以DSI、NDSI和RSI為自變量，建立冬小麥全蝕病病情指數高光譜監(jiān)測模型。通過建模得出3個模型的決定系數R2均在0.65以上，其中DSI(R818,R534)與病害指數的相關性最高，R2=0.860 5，RMSE=7.3%，MAE=19.1%；其次為RSI(R530,R758)，R2=0.659 0，RMSE=11.5%，MAE=20.6%；NDSI(R758,R530)與病害指數的相關性最低，R2=0.653 6，RMSE=11.6%，MAE=38.9%。綜上可知，DSI(R818,R534)對冬小麥全蝕病病情指數的監(jiān)測精度最高，并采用同期獨立驗證數據對3個模型的監(jiān)測精度進行檢驗，驗證精度R2均在0.6以上，其中DSI(R818,R534)構建的模型驗證精度最高(R2=0.76，RMSE=14.9%，MAE=11.7%，驗證樣本數n=20)，結果如圖7所示。

2.2.3冬小麥全蝕病估測的偏最小二乘模型

為探討多元線性回歸模型對冬小麥全蝕病DI的估測精度，本研究利用偏最小二乘回歸分析(PLSR)的方法，選取與光譜指數法相同的建模樣本，通過3種光譜指數與病情指數進行模型的構建，得到多元線性回歸模型y=0.576 6-2.604 6DSI(R818,R534)-0.015 0NDSI(R758,R530)+0.390 3RSI(R530,R758)，模型預測值與實測值的散點圖如圖8所示(R2=0.629 2,RMSE=10.2%,MAE=16.6%)。由偏最小二乘回歸分析結果可得出，相較于DSI(R818,R534)差值光譜指數，NDSI (R758,R530)和RSI(R530,R758)對回歸模型的貢獻度均較低，DSI(R818,R534)對模型貢獻度最高，與2.2.1節(jié)光譜指數和病情指數之間相關性結果一致。所以，本文最終利用DSI光譜指數單獨與病情指數構建線性回歸模型y=-6.490 1DSI(R818,R534)+1.461 3(R2=0.860 5,RMSE=7.3%,MAE=19.1%)。

圖7 病情指數反演與驗證Fig.7 Retrieving and validation of disease index

圖8 利用PLSR得到的實測值與預測值的相關性Fig.8 Correlation between predicted and measured disease indexes based on PLSR

2.3 冬小麥全蝕病病情指數反演與精度驗證

經過校正拼接處理后的無人機高光譜影像見圖9。從影像中提取40個樣點的DSI(R818,R534) 以構建冬小麥全蝕病病害指數估測模型，并應用于灌漿期的無人機高光譜影像，依據病害分級標準，通過ArcGIS 10.2軟件對冬小麥全蝕病病情指數閾值進行劃分，得到研究區(qū)域冬小麥全蝕病病情指數等級空間分布圖，見圖 10。通過獨立驗證數據集對該模型進行驗證，得到模型預測值與實測值的驗證精度(R2=0.76，RMSE=14.9%，MAE=11.7%，n=20)。

圖9 高光譜影像及采樣點分布Fig.9 Hyperspectral imagery and samples distribution in study area

圖10 冬小麥全蝕病病情指數等級空間分布圖Fig.10 Spatial distribution of disease index of wheat take-all in study area

綜上可知，通過對PLSR、線性回歸模型的建模精度進行檢驗，可得驗證精度R2均在0.62以上，其中DSI(R818,R534)構建線性回歸模型的精度更高。因此本文采用DSI(R818,R534)與DI構建的線性回歸模型y=-6.490 1x+1.461 3進行研究區(qū)冬小麥全蝕病的病害等級反演。

3 結論

(1)在462～874 nm波長范圍內(4～102波段)，UHD185成像光譜儀與ASD手持地物光譜儀所獲取的不同病害程度的冬小麥全蝕病光譜特征具有良好的一致性，二者光譜數據相關性顯著，決定系數R2達0.97以上，說明其數據具有較高的精度與準確性。

(2)冬小麥全蝕病病情指數與冠層反射率在734～874 nm波長范圍內呈顯著負相關。病害指數越高，冬小麥冠層光譜反射率越低。

(3)將DSI(R818,R534)、NDSI(R758,R530)和RSI(R530,R758)作為自變量與病情指數構建偏最小二乘回歸模型，盡管參與估測的光譜指數更為豐富，但模型驗證結果為R2=0.629 2,RMSE=10.2%,MAE=16.6%，其中DSI(R818,R534)對模型貢獻度最高，然而對DI 的估測較光譜指數DSI(R818,R534)模型精度差，只能粗略估測病情指數。

(4)DSI(R818,R534)差異光譜指數構建的線性回歸模型為病情指數的最優(yōu)估測模型，模型精度為R2=0.860 5,RMSE=7.3%,MAE=19.1%，模型驗證精度為R2=0.76，RMSE=14.9%，MAE=11.7%。將構建的線性模型應用于無人機成像高光譜影像，并把反演后的病情指數分成4個等級，生成研究區(qū)冬小麥全蝕病病情指數等級空間分布圖，為來年防治冬小麥全蝕病精準拌種提供了底圖，也為冬小麥及其他作物病害在田塊尺度上的監(jiān)測提供了解決思路。