基于高光譜數據與針葉電阻率的松材線蟲病早期診斷研究

劉 飛,談家金*,謝谷艾,蘇 峻,葉建仁

1. 南京林業大學林學院,南方現代林業協同創新中心,江蘇 南京 210037 2. 江西省林業科學院林業有害生物防治研究所,江西 南昌 330013 3. 南京林業大學理學院,江蘇 南京 210037

引 言

松材線蟲(Bursaphelenchusxylophilus)是目前對我國森林危害最嚴重的林業檢疫性有害生物,由其引起的松材線蟲病是一種松樹毀滅性病害[1-3]。我國自1982年發現該病以來,病情經過41年不斷發展,目前全國共有19個省(自治區、直轄市)701個縣級行政區發生疫情,松材線蟲病的準確診斷對該病的有效防控至關重要。該病現行的診斷技術包括林間癥狀診斷法、病原線蟲鑒定法、流膠法等,這些方法受到條件或技術的限制而不完善;為了能夠有效發現和控制松材線蟲病,其爆發前的早期診斷更是成為當前研究的熱點[4]。如高瑞賀等發現健康馬尾松感染松材線蟲病后,其針葉的蒸騰速率和凈光合速率變化顯著,可作為松材線蟲病早期診斷的指標[5];Gaspar等在接種松材線蟲后,發現松樹3-蒈烯揮發量相比健康松樹顯著提高,可用于松材線蟲病早期診斷的指標[6]。

大量研究表明,植物在受生物脅迫或非生物脅迫過程中,其自身色素、含水量等生理生化指標會發生一定的變化[7-8]。而植物對電磁輻射的吸收與反射特性會隨著這些生理生化參數的改變而改變,在遙感上表現為不同的光譜反射率。因而遙感在研究植物脅迫方面應用廣泛,可見光區和近紅外區是指示植物色素含量變化和健康情況的關鍵波段[9-11]。健康馬尾松在感染松材線蟲病后,其葉片內部結構會被破壞,光合色素和水分含量也會發生變化[12];伴隨病害不斷的發展,馬尾松的光譜反射率呈現出較大變化[13]。松材線蟲病在病情發展的早期階段,感病松樹自身不會發生明顯癥狀,且葉綠素和水分含量變化不明顯,但當肉眼可以觀察到針葉葉色明顯變色時,感病松樹會在短期內萎蔫死亡。目前研究中,無論是基于光譜本身的直接診斷方法,或是基于生理生化組分的間接診斷方法在病害的早期階段效果并不理想[14-15]。

為了實現對感病植株的早期診斷,可以考慮尋找一種基于光譜的新方法用于判斷植株是否感病,這將對松材線蟲病的防治具有重要的實踐意義與應用價值。馬尾松感染松材線蟲病后除生理生化指標發生變化外,其生物電也必然會發生一定的變化[16]。本研究擬從電學角度出發,在研究感病松樹針葉電阻率變化的基礎上,將基于馬尾松針葉反射光譜數據建立針葉電阻率的預測模型,以實現在感病后的早期階段對馬尾松松材線蟲病的快速判別,為利用無人機、衛星遙感對松林進行松材線蟲病監測提供理論依據。

1 實驗部分

1.1 松材線蟲接種、病情觀察及樣品采集

實驗地設在安徽省休寧縣商山鎮五方水庫的一處馬尾松林內(118°21′63″E,29°69′23″N;156.7 m)。2021年8月17號采用注干法對8～9年生馬尾松接種松材線蟲,實驗共設置處理和對照兩組,處理組和對照組的馬尾松分別設置6棵和4棵。在松樹主干距離地面1.5 m處,使用直徑6 mm的鉆頭(鉆頭與主干夾角約45°)鉆入木質部4 cm,在鉆孔內放入適量的無菌棉球;使用1 mL的槍頭(槍頭前端已被修剪為注射器針頭狀,增加接觸面積)插入鉆孔,用微量進液器注入線蟲懸浮液,處理組接種線蟲量為30 000條·株-1(濃度30 000條·mL-1),對照組接等量無菌水。接種后每4～7天對馬尾松發病狀況進行觀察并記錄感病指數。馬尾松病情分級與感病指數按談家金等的方法[17],略加修改后計算:0級正常,針葉綠色,代表值為0;1級為1/2以下針葉褪綠、1/4以下葉發黃,代表值為1;2級為1/2以上葉褪綠、1/4～3/4葉發黃,代表值為2;3級為3/4以上葉發黃、1/2以下葉變紅,代表值為3;4級為1/2以上葉變紅、植株瀕死或死亡,代表值為4。計算公式:感病指數=Σ(各病級株數×病級代表值)×100/(總株數×最高病級的代表值)。接種當天以及接種后7、11、15、19和26 d,采集外觀表現健康的針葉進行光譜和電阻的測定。樣品采集時間設置在上午10:00—12:00,選擇晴朗無風無云或少云天氣,從馬尾松向陽面3個位置(即冠層上、中、下3個位置各取一個點)采集當年生成熟針葉,將采集好的樣本針葉用密封袋封口,放置在冰盒內帶回室內進行光譜和電阻的測定。

1.2 針葉光譜反射率測定

采用Ocean Optics USB2000+光譜儀對針葉光譜反射率進行測定,該光譜儀由美國Ocean Insight公司生產,其作用波段范圍為350～1 000 nm,光譜分辨率為1.5 nm,光譜采樣間隔為0.44 nm。光譜儀探頭視場角25.4°,探頭距離樣本表面距離8 cm,采用10 W鹵素燈為光源,方位角為45°(光源與樣本表面的夾角)。每次測定前,使用白色面板進行光譜校正。測定儀器各配件如圖1所示,開關電源(a)、探頭(b)、樣本槽(c)和光源(d)均放置在內部密閉黑色的箱子內,如圖1(A)所示。測定過程中設定光譜平均次數為10次,積分時間100 ms。對針葉樣本進行測定時,使用黑色雙面膠將針葉中間部分整齊平鋪在黑色樣本板上,并剪去多余部分,針葉之間不能有縫隙,如圖1(B)所示。將冠層上、中、下3個位置光譜反射率平均值作為該植株的光譜反射率。

圖1 針葉光譜反射率測定儀器各配件(A)黑色箱子內部各組件:(a)開關電源;(b)探頭;(c)樣本槽;(d)光源;(B)樣品板與待測針葉樣本Fig.1 Accessories of needle spectral reflectance measurement instrument(A) Internal components of the black box:(a) Switching power supply;(b) The probe;(c) The sample slot;(d) The light source;(B) Sample plates and sample of needles to be tested

1.3 針葉電阻率測定

采用M4070 LCR測試儀測定針葉電阻,該儀器是由中國晶研儀器科技有限公司生產,其量程范圍為0.00 Ω～10.000 MΩ,基本測量準確度0.5%。使用3D打印機打印黑色薄板,將兩個電極插入對應位置,在半球形凹槽內放入導電膏(濟寧高新區金諾特醫用凝膠廠生產),方便針葉電阻的測量,如圖2所示。將針葉橫截面近似看成半橢圓形,剪取針葉中部4 cm,使用游標卡尺(分辨率0.01 mm)對針葉的寬度和厚度進行測量。將針葉與M4070 LCR測試儀連接后3～5 s的讀數作為其電阻值,從而計算出其電阻率。冠層向陽面上、中、下3個位置各取5根針葉進行電阻的測定,將5根針葉電阻率平均值(電阻率異常值的樣本除外)作為對應冠層位置針葉電阻率,將上、中、下3個位置電阻率平均值作為該植株的針葉電阻率。

圖2 針葉電阻測定過程Fig.2 Needle resistance measurement process

1.4 數據分析

1.4.1 數據預處理

地面光譜與針葉電阻同步采樣共6次,將獲取的數據去除明顯異常值后,共獲得52組數據。其中,處理組含34組數據,對照組含18組數據。使用Excel數據分析工具對光譜數據和電阻數據處理匯總;使用SPSS軟件對電阻率數據進行分析,并在Origin 2018中繪圖。

在MATLAB R2021b中調用Savitzky-Golay函數對光譜數據進行平滑處理,既可以消除數據中的噪聲部分,又能將數據中重要的信息保留,提高了數據的平滑性。將平滑后的光譜定義成原始光譜(OR)。光譜數據自身具有復雜性和高維數特點,為了提取更多有用的信息,需要對光譜數據進行一些復雜的數學變換,為預測模型的創建提供充足的信息源。目前通常的光譜變換包括:一階導數(FD)、二階導數(SD)、對數變換(LOG)、倒數變換(1/R)和連續統去除法(CR)。

1.4.2 特征波段選擇

高光譜數據具有高維、復雜的特點,為了篩選出與針葉電阻率相關性高的波段,采用隨機森林進行特征波段篩選。隨機森林(random forest,RF)是一種非線性模型,以決策樹作為基學習器,使用裝袋(Bagging)原理來處理訓練數據集[18]。由于其計算過程相對簡單,對樣本和變量的選擇過程是隨機的,因而能夠在多種領域應用并且均可以得到很好的分類效果[19]。這些特點使得RF特別適用于光譜數據的研究。本次實驗使用RF進行特征波段的篩選,RF中樹的數量nTrees設置300,輸入的特征數量10。

1.4.3 最小二乘支持向量機回歸模型

最小二乘支持向量機(least squares support vector machine,LSSVM)算法是基于支持向量機算法的一種改進,極大提高了求解問題的速度和收斂精度[20],且不再受到樣本數量的影響,具有很強的泛化能力。

采用LSSVM算法對針葉光譜反射率與針葉電阻率進行回歸分析,從而對馬尾松針葉電阻率進行反演。以光譜變換下的特征波段為自變量,針葉電阻率為因變量,基于LSSVM算法構建馬尾松針葉電阻率預測模型。為了使得模型測試的結果更客觀,80%的樣本選作建模集,20%的樣本用作驗證集。主要用以下三個指標評估馬尾松針葉電阻率預測模型的精度:決定系數(coefficient of determination,R2)、均方根誤差(root mean square error,RMSE)、平均絕對誤差(mean absolute error,MAE)。當預測模型的R2趨近1、RMSE和MAE越小時,表明其精度越高。

2 結果與討論

2.1 接種后馬尾松病情發展

接種后11 d,發現有1棵接種線蟲的馬尾松1個2級側枝針葉開始褪綠,此時,病情指數為4;接種后15 d,該側枝針葉變為黃褐色,此時,病情指數為4;接種后19 d,該側枝針葉變成紅褐色,針葉已枯死,此時,已有3棵接種線蟲的馬尾松針葉變色,病情指數為12.5。接種后26 d,已有4棵接種線蟲的馬尾松極少數針葉變色,病情指數為16.7。接種后67 d,接種線蟲的6棵馬尾松少數針葉均已變色,病情指數為25。

2.2 接種后馬尾松針葉電阻率的變化

接種實驗結果顯示,馬尾松顯癥早期,針葉電阻率顯著升高。由圖3可知,接種后15 d,此時僅1棵接蟲馬尾松的1個側枝針葉變色,接蟲株和對照株的電阻率開始達到極顯著差異,這種差異一直維持到接種后26 d(p<0.01)。整個實驗期間,相比接蟲株的針葉電阻率,對照株的針葉電阻率波動較小。因此,從電阻率的大小可以初步判斷當馬尾松針葉電阻率大于323 Ω·m時,該松樹已經感病。

圖3 接種松材線蟲后馬尾松針葉電阻率的變化**表示在0.01水平差異顯著Fig.3 Changes in needle resistivity of Pinus massoniana inoculated with B. xylophilus** Shows significant correlation difference at 0.01 level

2.3 接種后馬尾松針葉光譜數據特征分析

由于部分波段區域出現大量噪聲,主要對490～980 nm范圍內的馬尾松針葉反射光譜進行分析。為了豐富基于光譜建立針葉電阻率模型的信息,使用5種數學變換方式對所有有效光譜數據進行光譜預處理,僅對數變換沒有改變光譜曲線的趨勢[圖4(a—f)]。

圖4 不同預處理方法的光譜變化(a):OR;(b):FD;(c):SD;(d):LOG;(e):1/R;(f):CRFig.4 Spectral changes of different pretreatment methods(a):OR;(b):FD;(c):SD;(d):LOG;(e):1/R;(f):CR

由表1可知,通過隨機森林算法對各預處理方法進行特征波段篩選后發現存在較大不同。原始光譜(OR)的特征波段集中在近紅外區域,一階導數和連續統去除法光譜的特征波段集中分布在綠光和近紅外區域,對數變換和倒數變換篩選出的波段主要分布在紅光和近紅外區域。綜合比較發現,特征波段主要分布在近紅外區域,說明病害發生后的早期針葉內細胞結構可能發生了變化,這與高瑞賀等研究的感染松材線蟲病引起馬尾松針葉光合作用相關機構改變的結果相一致[5]。大量研究表明,在可見光和近紅外區域的反射光譜在病蟲害監測中具有很高的應用價值,它們通常與病害的發展過程密切相關[9-10]。

表1 基于隨機森林算法篩選的特征波段Table 1 Feature bands screening based on random forest algorithm

2.4 馬尾松針葉電阻率預測模型分析

基于原始光譜進行五種光譜變換,使用隨機森林對原始光譜和各預處理后的光譜數據進行篩選,從而得到與針葉電阻率相關特征波段,最終利用LSSVM算法建立的針葉電阻率預測模型。該模型的預測效果匯總于表2,發現不同預處理間預測效果的綜合表現具有差異,其中二階導數模型預測效果綜合表現最好,建模集與驗證集的平均R2、RMSE和MAE分別為0.848、32.332和7.067。

表2 基于不同光譜變換篩選輸入變量的模型表現Table 2 Performance of model based on the different spectral transformation input variable selection

各預測模型間除了倒數變換和連續統去除法變換外,其他光譜變換后的預測模型比原始光譜預測模型驗證集的R2均有所提升,RMS和MAE均有所降低。一階導數預測模型的改善效果較小,建模集與驗證集的平均R2僅提升1.5%,但RMSE和MAE卻分別升高28%和1.7%。而二階導數模型預測效果提升較為明顯,R2提升了4.2%,RMSE和MAE分別降低2.5%和18.9%。對數模型預測效果也有所提升,R2提升了3.1%,但RMSE和MAE分別升高了57.3%和11.3%。綜合比較,基于二階導數選擇的特征變量在最小二乘支持向量機算法預測模型中綜合表現最優,特征波段為594.986、646.107、646.451、782.896、784.841、839.164、863.890、902.021、947.901、962.315;其次為LOG、FD、1/R和CR,這表明對光譜預處理的必要性,處理后的光譜數據增強了樣本的光譜反應特征,進一步提升反演模型的精度。

圖5為基于二階導數(SD)變換的LSSVM建模方法反演結果散點圖,從圖5中可以清晰看出,該模型的實測值樣點和預測值樣點多數分布在y=x線附近,能夠較好地對馬尾松針葉電阻率進行估測。因此,認為SD-RF-LSSVM可作為馬尾松針葉電阻率大小的最優估測模型。相較于傳統的林間癥狀診斷法、病原線蟲鑒定法以及流膠法[4],本文提出的基于SD-RF-LSSVM模型的光譜檢測方法不僅可以實現對該病害的快速無損診斷;此外,該方法能夠在病樹極少數針葉變色時就完成對病害的監測,相比之下,傳統的遙感監測方法需要等到病樹針葉大部分已經明顯變色才能發現病害[16]。因此,該方法可以在病樹出現極少數針葉變色時進行早期預警和監測,及時發現病害的發生和變化情況,從而更好地制定病害防治措施。

圖5 基于二階導數數據變換的最小二乘支持向量機模型檢驗結果Fig.5 Test results of LSSVM model based on the SD data transformation

3 結 論

多年來,我國在松材線蟲病的致病機理及防治方法上取得許多重要進展,然而在該病的早期診斷方面一直未有突破,而利用遙感技術對該病的早期監測和預警長期處在探索階段。本研究表明,使用針葉反射光譜反演針葉電阻率是可行的,且SD-RF-LSSVM建立的預測模型精度最高,可用于針葉電阻率的快速估測,為實現基于遙感的松材線蟲病早期診斷與監測提供了思路與方法。

實驗僅對針葉葉片層面的光譜特征進行分析,雖可為將來冠層尺度研究提供參考,但冠層尺度的松樹冠層監測更復雜,光譜識別與分析受到冠層形狀、土壤、植被等因素影響。需要在冠層尺度比較健康和感病馬尾松的光譜特征差異,從而為使用遙感技術實現病害的快速有效識別提供方法。此外樣本數量較少,模型需要進一步檢驗和驗證。同時,環境條件對結果的影響需要進一步研究,增強模型的穩定性和適用性是今后研究的重點。