基于Faster R-CNN的松材線蟲病受害木識別與定位

徐信羅 陶 歡 李存軍 程 成 郭 杭 周靜平

(1.南昌大學信息工程學院， 南昌 330031； 2.北京農業信息技術研究中心， 北京 100097)

0 引言

松材線蟲病是松樹的毀滅性病害，具有很強的傳染性，主要通過傳播媒介松墨天牛[1](Monochamusalternatus)進行擴散。松樹一旦感染此病(以下稱為“受害木”)，從出現癥狀到死亡只需40 d左右，整片松林的毀滅只需3～5年[2]。KIYOHARA等[3]進行了大規模的接種試驗，于1971年首次證實了松材線蟲(Pine wilt nematode，PWN)對松樹的致病性。在分子植物病理學中，松材線蟲是十大植物寄生線蟲之一[4]。松樹染病后，其針葉由有光澤的綠色變為黃褐色乃至紅褐色，直至最后枯萎死亡[5]。松材線蟲病起源于北美洲，之后傳播至日本[6]。我國自1982年在南京中山陵首次發現感染松材線蟲病的黑松[7]后，該病在我國不斷擴散，并造成大量松樹死亡。據國家林業和草原局2019年2月1日發布的松材線蟲病疫區公告，全國范圍內松材線蟲病疫區已擴散至18個省的588個縣，新增疫區282個，發生面積64.93萬hm2[8]。據不完全統計，30多年來，全國因松材線蟲病致死的松樹累計達數十億株，造成了上千億元直接和間接的經濟損失。

殺滅松材線蟲病的傳播媒介松墨天牛是目前防控松材線蟲病的主要手段，最有效的措施是直接對受害木進行砍伐，然后對其進行熏蒸、焚燒或微波加熱等無害化處理[9]，也可采用懸掛誘捕器[10-11]和投放天敵的方式。疫情調查是松材線蟲病防控的工作基礎，現有受害木的監測手段主要有地面調查[12]、衛星遙感監測[13-14]和無人機遙感監測[15-20]。地面調查和衛星遙感監測有其局限性。無人機遙感具有靈活性高、應用周期短、時間和空間分辨率高、成本低、操作簡便等優點[21]，目前已在森林資源調查、森林火災監測、森林病蟲害監測防治、森林信息提取等方面得到廣泛應用[21-23]。在受害木的檢測中，使用無人機遙感既能節省大量的人力和物力，又能克服空間分辨率的限制，達到監測單株受害木的目的。

基于無人機的受害木識別已有大量的研究報道。李衛正等[15]用無人機獲取高空間分辨率的松林影像，正射處理后導入Geo Link軟件，通過目視判讀尋找受害木。呂曉君等[16]根據感病松樹樹冠顏色的變化，對無人機采集的數字正射圖像進行目視判讀。利用人工對無人機影像中的受害木進行判讀，效率低且主觀性強。陶歡等[17]對獲取的無人機影像采用HSV(色調、飽和度、明度)閾值法實現變色松樹的識別，能有效提高人工判斷的效率。隨著圖像分析技術的發展，劉遐齡等[18]用無人機獲取高分辨率影像，采用多模板識別法對不同染病階段的受害木進行識別，結果表明，相比于目視判讀，模版匹配方法能有效提高受害木的檢測效率。當機器學習成為研究熱點后，一些研究者嘗試使用機器學習的方法對受害木進行檢測[19-20]。

基于深度學習的Faster R-CNN[24]目標檢測算法在番茄[25]、船[26-27]、鳥[28]和飛機[27-29]等物體的檢測中都取得了較好的效果。Faster R-CNN是一種基于區域建議的目標檢測算法，在Fast R-CNN[30]和R-CNN[31]算法基礎上改進后提出。與基于回歸的SSD[32]和YOLO[33]等目標檢測算法相比，Faster R-CNN算法的檢測精度更高[34]。目前，在受害木的識別中尚未見使用深度學習的方法。本文提出基于深度學習的高空間分辨率無人機遙感影像的受害木自動檢測方法，同時考慮因其他原因(干旱或自然死亡等)致死的松樹(以下稱為“其他枯死樹”)和紅色闊葉樹對受害木識別的影響，使用包含受害木、其他枯死樹和紅色闊葉樹的數據集對模型進行訓練，以達到提高受害木識別效果和檢測精度的目的。

1 研究區概況

研究區位于福建省晉江市紫帽鎮(24.87°～24.91°N，118.47°～118.53°E)，如圖1a所示，該地區為松材線蟲病重疫區[8]。所選區域的總面積為4.25 km2，最高海拔為517.8 m，如圖1b所示。晉江市屬于亞熱帶季風氣候，夏季長冬季短，年平均光照時間為2 100 h，年平均溫度為20～21℃，七、八月是溫度最高的2個月(平均氣溫27.8～29.4℃)。年平均降水量為911～1 231 mm，降雨主要集中在夏季，夏季以西南風為主，其他季節盛行東北風。在現有森林資源中，松林面積16 km2，占全市林分總面積的26.15%。晉江市于2015年秋季和2016年春季開展松材線蟲病及受害木普查，紫帽鎮為重點調查區域且該區域的疫情較為嚴重。由于寄主植物(松樹)和寄主昆蟲(松墨天牛)的存在以及舒適的生存環境，松材線蟲病在該區域呈爆發式傳播。

圖1 研究區域Fig.1 Research area

2 材料與方法

2.1 數據獲取與預處理

使用QK-4071型無人機(配備包含紅、綠、藍3個光譜的SONY A5100型數碼相機)獲取研究區域的影像，將全球導航衛星系統(Global navigation satellite system, GNSS)和慣性測量單元(Inertial measurement unit, IMU)模塊整合到無人機平臺，確保水平和垂直方向的位置誤差在2 m和5 m左右。數碼相機的鏡頭與位置和方向系統(Position and orientation system, POS)信息的記錄同步，以確保將POS信息附加到每幅圖像上。

首先從Google Earth衛星影像中了解研究區域的地理情況并根據天氣情況制定相應的飛行計劃。然后在2018年3月11—13日利用無人機獲取影像，在晴朗無風條件下共飛行了13架次，分別為：3月11日11:00—16:00，飛行5架次；3月12日11:00—16:00，飛行4架次；3月13日10:00—13:00，飛行4架次。考慮到實際地形、植被狀況以及所需覆蓋的區域，將無人機飛行的相對高度設置為300 m，每次飛行的航向重疊率和旁向重疊率分別為60%和50%。每次飛行后，檢查每幅影像的數據質量，將無效影像刪除并再次拍攝。最后使用Pix- 4D軟件進行圖像校準和鑲嵌，得到空間分辨率為8.47 cm的數字正射模型(Digital orthogonal model, DOM)影像，利用ENVI軟件對影像進行增強。使用無人機自身攜帶的差分自動駕駛儀可以精確給出每幅影像中心點經緯度坐標，再結合無人機航向、姿態數據和相機參數，可精確給出影像上每個像元的坐標信息，為受害木的定位提供準確位置信息。

2.2 數據集準備

受害木表現出紅褐色的顏色特征，在無人機影像上與綠色的健康松樹存在明顯的色差。通過2018年4月1—3日的地面調查以及對無人機影像的分析發現，研究區域存在一些顏色、紋理與受害木相近的樹木，比如其他枯死樹，此類松樹以灰褐色為主；紅色闊葉樹呈黃褐色，此類樹木與枯死樹顏色最為相近，3種樹木的無人機影像和實地圖像如表1所示。其他枯死樹和紅色闊葉樹沒有松材線蟲病，故無需對其進行處理，但這兩種樹木對受害木的識別會造成一定干擾。

表1 3種樹木的無人機影像和實地圖像Tab.1 UAV images and field photos of three trees

首先將研究區劃分為訓練區和測試區，如圖1b所示，藍色區域為訓練區，紅色區域為測試區。選擇較大的區域作為訓練區，是為了保證深度學習模型學習到足夠多的特征信息。然后對訓練區和測試區進行影像裁剪，人工對訓練區內的受害木進行標注并裁剪圖像(256×256，圖像尺寸單位為像素，下同)，在訓練區共采集到1 283幅圖像，按8∶2的比例將訓練區的圖像隨機分成訓練集和驗證集，訓練集包含1 026幅圖像，驗證集包含257幅圖像。通過ArcGIS軟件中的Split Raster工具將測試區的圖像裁剪為256×256，共采集到279幅圖像。由于其他枯死樹和紅色闊葉樹的存在，故將這兩種樹木與受害木同時用于模型訓練，則該數據集包括受害木、其他枯死樹和紅色闊葉樹3類樹木。在使用LabelImg軟件制作與圖像一一對應的可擴展標記語言(XML)文件時，將3種類型的樹木都做上標記，將該數據集稱為數據集1。此外，建立只包含受害木信息的數據集，即在制作XML文件時只對圖像中受害木進行標記，將該數據集稱為數據集2，用于驗證在考慮了其他枯死樹和紅色闊葉樹后是否有助于提高受害木的檢測精度。數據集1中各類樹木的統計情況如表2所示，除受害木外，研究區內的其他枯死樹數量相對較多，而紅色闊葉樹數量較少，整個區域共標記2 102棵受害木、778棵其他枯死樹和192棵紅色闊葉樹，而數據集2只包含表2中的受害木部分。

2.3 研究方法

2.3.1受害木識別算法的整體框架

Faster R-CNN的整體框架如圖2所示，其中FC表示全連接層，Conv表示卷積層。Faster R-CNN目標檢測算法主要分為兩部分，分別是區域生成網絡(Region proposal network，RPN)和Fast R-CNN。Faster R-CNN通過共享卷積的方式將這兩部分連接起來。使用VGG16[35]作為特征提取網絡并使用線性整流函數(Rectified linear unit，ReLU)作為激活函數，VGG16包含13個卷積層、5個最大池化層和3個全連接層。當輸入的無人機遙感圖像經過VGG16特征提取網絡后，會對受害木、其他枯死樹和紅色闊葉樹的特征信息進行提取，并輸出特征圖。RPN網絡以卷積神經網絡輸出的特征圖作為輸入，并通過滑窗(尺寸為3×3)獲得錨框，錨為每個滑窗的中心，結合不同尺寸和比例的區域建議，每個錨產生9個不同的錨框，然后輸出一系列可能包含受害木、其他枯死樹和紅色闊葉樹的矩形候選框和相應的候選框得分。由于RPN產生的區域候選框的尺寸不同，故使用RoI Pooling層將不同尺寸的區域候選框映射成固定尺寸后輸入全連接層，RoI Pooling以特征提取網絡輸出的特征圖和RPN網絡輸出的區域候選框作為輸入。最后利用Softmax層對每個候選框進行分類并輸出所屬類別的得分；同時利用邊框回歸獲得每個候選框相對實際位置的偏移量預測值，用于修正候選框的位置，以得到更精確的受害木邊界框。

圖2 Faster R-CNN整體框架Fig.2 Overall framework of Faster R-CNN

圖2 Faster R-CNN整體框架Fig.2 Overall framework of Faster R-CNN

在原始的RPN網絡中，錨框是基于多尺寸({1282,2562,5122})和多比例({1∶1, 1∶2, 2∶1})產生的。PASCAL VOC數據集中的圖像尺寸和目標尺寸都比較大，則原始的錨框尺寸適合于該數據集，并能取得較好的實驗結果。根據實際情況對數據集中受害木的冠幅面積進行統計分析，由于3類樹木的冠幅尺寸相近，故只對受害木的冠幅進行統計。由圖3可以看出，冠幅面積主要在102～602之間，且在302左右最為集中，原始錨框的像素面積尺寸遠遠超過了本數據集中受害木冠幅的尺寸，因此將錨框的尺寸修改為{162,322,642}，使之能達到更好的識別效果。

圖3 冠幅面積統計Fig.3 Crown area statistics

2.3.2坐標轉換與點位合并

為了方便森林防護人員及時了解受害木的分布情況并進行高效率的處理工作，需給出模型預測結果中受害木的平面坐標信息并以點的形式在測試區域中標出。模型的輸出結果包含預測邊界框的像素位置信息，即相對于圖像左上角的像素坐標信息，結果用4個數值表示預測邊界框的位置，分別是左上角和右下角的像素坐標值ximin、yimin、ximax、yimax。因此需要將預測邊界框中心點的像素坐標轉換為平面坐標點信息，然后用ArcGIS軟件進行標點，所選用的坐標系統為WGS_1984_UTM_zone_50N。

圖4a中的XY坐標系為平面坐標系，xy坐標系為像素坐標系，大的矩形框表示圖像，小的矩形框表示圖像中的某一預測邊界框。結合圖像獲取和圖像裁剪過程，可知每幅圖像中心點(X0,Y0)的平面坐標，然后根據(X0,Y0)計算每個預測邊界框的中心點的平面坐標。首先計算圖像左上角的平面坐標值

圖4 坐標轉換與點位合并Fig.4 Coordinate transformation and point combination

(1)

其中

L=b=21.69 m

式中L——圖像長度所表示的實際地理長度，m

b——圖像寬度所表示的實際地理寬度，m

(X01,Y01)——圖像左上角的平面坐標,m

然后計算圖像中預測邊界框中心點像素坐標值

(2)

式中i——圖像上預測框的序號

(xi,yi)——邊界框中心點像素坐標

最后將預測框中心點的像素坐標轉換為平面坐標

(3)

其中

Δx=Δy=8.47 cm

式中 Δx——水平方向上的空間分辨率，cm

Δy——垂直方向上的空間分辨率，cm

(Xi,Yi)——第i個預測邊界框中心點平面坐標，cm

由此可得到每個預測邊界框中心點的平面坐標。

在對無人機影像進行裁剪的過程中可能將某一棵受害木切割到兩幅甚至更多幅子圖像中，故在定位結果中會存在某一棵受害木同時由兩個以上點位表示的情況。為了避免這種情況，需要對同時表示一棵受害木的多個點進行合并。通過對整個測試區域的分析發現，受害木的分布較為離散，不存在連片的受害木；對受害木的冠幅尺寸進行分析，由圖4b可以看出，冠幅的長、寬主要集中在80像素以內，再結合空間分辨率8.47 cm，故將點位合并的歐氏距離閾值設置為3.40 m。合并后的點位坐標值為原來各點位的平均坐標值。

2.4 Faster R-CNN模型訓練

2.4.1實驗平臺和參數設置

實驗操作平臺為Ubuntu16.04計算機，Intel Core i7-3770 CPU@3.40 GHz，8位英特爾處理器，NVIDIA GeForce GTX 1080TiGPU，使用TensorFlow作為深度學習框架。在參數設置方面，迭代次數設置為40 000次，初始學習率設為0.001，每次迭代訓練圖像的數量(batch size)為256，學習率的衰減系數(γ)和動量分別為0.1和0.9[25]，非極大值抑制(NMS)的閾值為0.7[29]，置信度閾值設置為0.8。

2.4.2評價指標

為了定量評價模型的性能，采用交并比(Intersection over Union，IoU)、準確率P(Precision)、召回率R(Recall)、總體精度F1值和平均精度(Average precision)作為評價指標。

(1)交并比

交并比定義為模型預測邊界框與真實標記框的重疊比率，用于評價目標是否被準確預測，結果越趨近于1表示預測越準確。計算公式為

(4)

式中I——交并比A——面積函數

Bpred——預測邊界框的區域

Btruth——真實標記框的區域

(2)召回率和準確率

R和P的計算依賴于正確檢測(TP)、錯檢(FP)和漏檢(FN) 3個參數。如果某一預測邊界框與同類的真實標記框的交并比大于某一確定的閾值(本研究中I設定為0.5)[27]，則將該預測邊界框定義為TP，即為正確檢測。否則將該預測框定義為FP，即為錯檢。如果某個真實標記框找不到相應的預測邊界框與之相匹配，則將這一真實標記框定義為FN，即為漏檢。P和R的計算公式為

(5)

(6)

式中TP——正確檢測的樣本數量

FP——錯誤檢測的樣本數量

FN——漏檢樣本數量

P表示在某個類別的所有預測結果中，正確檢測的比例。R表示在某個類別的所有真實標記框中，正確檢測的比例。

(3)總體精度和平均精度

F1值[35]用于評價模型的整體性能，計算公式為

(7)

式中F1——總體精度

F1值為P和R的調和平均數，在0(最差)和1(最佳)之間變化。

平均精度是評價模型性能的重要指標，計算公式為

(8)

式中IAP——平均精度

若一個模型在不同的R下都能保持較高的P，則IAP就越高，說明模型對此類檢測的表現較好。

3 實驗結果

3.1 受害木識別

將由數據集1訓練得到的模型對測試集進行測試，得到的受害木IAP為78.42%，其他枯死樹和紅色闊葉樹的IAP分別為63.27%和33.20%，紅色闊葉樹的IAP較低，可能與其樣本量相對較少有關。本研究以受害木的識別為重點，由于其他枯死樹和紅色闊葉樹不帶有松材線蟲病，且不需要森林防護人員進行及時清理，故對這兩種樹木的檢測結果不做具體分析。

該模型檢測的各類樹木的準確率、召回率和F1值的計算結果如表3所示，模型從572棵受害木中正確檢測出504棵，漏檢了68棵，其召回率達到了88.11%，受害木漏檢率為11.89%；模型誤判了147棵受害木，原因是將其他枯死樹和紅色闊葉樹誤判為受害木，其準確率為77.42%；模型識別受害木的總體精度為82.42%。

表3 準確率、召回率和F1值的計算結果Tab.3 Calculation results for precision, recall and F1 scores

模型訓練過程中損失值(Loss)隨迭代次數變化的曲線如圖5所示。模型的迭代次數為40 000次，在模型訓練的初始階段其Loss值保持在一個較高的水平，隨著迭代次數的增加，Loss值存在一個緩慢振蕩下降的過程，當迭代到25 000次之后，Loss值基本趨于穩定且不再下降。總體來說，Loss曲線相對平滑且最終Loss值保持在0.14左右。

圖5 Faster R-CNN訓練時的損失值曲線Fig.5 Loss curve of Faster R-CNN training

3.2 受害木定位

從模型的預測結果中篩選出所有被正確檢測的受害木共504棵，然后將每個受害木預測邊界框的位置信息通過式(1)～(3)轉換得到其在平面坐標系中的點位坐標，并將點位坐標通過ArcGIS軟件在測試區域以點的形式標出。最后，通過設定的點位合并閾值將表示同一棵受害木的多個點位進行合并，最終正確定位出494棵受害木，點位合并過程如圖6a所示。最終的受害木定位結果如圖6b所示。

圖6 受害木定位合并過程與定位結果分布圖Fig.6 DPT positioning process and result distribution map

4 討論

使用只包含受害木信息的數據集2進行模型訓練與檢測。此外，使用原始錨框尺寸({1282,2562,5122})對2個數據集進行檢測實驗，所有實驗使用的平臺和參數設置都一致。表4為在不同模型下測試得到的IAP，表5為不同模型下受害木檢測的準確率、召回率和F1值。對比表4中模型1、3和模型2、4發現，在相同的錨框尺寸下，數據集1訓練模型比數據集2訓練模型的受害木的IAP分別高5.38個百分點和6.14個百分點；同理，表5中模型1、2比模型3、4的F1值分別提高了5.22個百分點和6.10個百分點，此外準確率和召回率都有不同程度的提高；結果表明，由數據集1訓練的模型在受害木識別上有更好的表現和檢測精度。對比表4中模型1、2和模型3、4發現，在相同的數據集下，根據實際目標大小將錨框尺寸修改為{162,322,642}后，受害木的IAP分別提高了2.68個百分點和3.34個百分點，且模型1中其他枯死樹和紅色闊葉樹的IAP均高于模型2，該結果與REN等[29]在小目標檢測中根據真實目標大小修改錨框尺寸后得到的實驗結果一致；同理，表5中模型1、3比模型2、4的受害木的F1值分別提高了0.68個百分點和1.56個百分點，準確率和召回率同樣有小幅度的提高。結果表明，在使用基于深度學習的目標檢測技術對受害木進行檢測時，根據受害木冠幅大小修改錨框尺寸后并將松林中存在的其他枯死樹和紅色闊葉樹的樣本信息加入數據集并用于模型訓練，能提高受害木識別效果和檢測精度。

表4 不同模型的IAPTab.4 IAP values of different models %

表5 受害木在不同模型下的準確率、召回率和F1值Tab.5 Precision, recall and F1 score values for DPT under different models %

使用目標檢測技術對受害木進行識別能達到較好的檢測精度，同時也會引起漏檢和誤判。最理想的檢測結果如圖7a所示，圖中5棵受害木均能被正確檢測出。圖7b中有2棵受害木，模型正確檢測出了1棵，漏檢了1棵。圖7c顯示的是將其他枯死樹誤判為受害木的情況，除了將1棵其他枯死樹誤判為受害木外，模型正確檢測出了其余受害木。圖7d是將紅色闊葉樹誤判為受害木的情況，模型將圖中受害木和其中一棵紅色闊葉樹正確檢測出來的同時，將另一顆紅色闊葉樹誤判為受害木。誤判的原因可能是其他枯死樹和紅色闊葉樹的顏色特征相近。圖7中deadpinetrees表示受害木，otherdeadtrees表示其他枯死樹，redbroadleavedtrees表示紅色闊葉樹，黃色矩形框表示漏檢，黑色矩形框表示誤判。

圖7 模型1檢測結果可視化Fig.7 Model 1 detection results visualization

使用深度學習的方法對受害木進行識別，F1值能達到82.42%，能夠滿足森林防護人員對受害木實地砍伐的要求。與傳統的目視判讀方法相比該方法能極大地提高受害木的識別效率并且能做到自動識別。與陶歡等[17]提出的HSV閾值法的總體精度(58%～65%)相比，檢測精度有明顯的提高。本研究的精度略低于劉遐齡等[18]提出的模板匹配方法的正確率(83.9%)，但本研究考慮了其他枯死樹和紅色闊葉樹的存在，精度的可靠性更高。使用深度學習方法檢測受害木能充分發揮計算機識別速度快的優點，從而達到及時檢測的目的。最后通過坐標轉換的方法對預測結果中的受害木進行定位，能提供每棵受害木精確的位置信息及總體分布情況，為森林防護人員尋找和清除染病松樹提供了技術支持，且能有效提高工作效率。本研究利用無人機遙感技術對松材線蟲病進行監測并使用目標檢測算法對染病松樹進行識別，表明通過結合無人機遙感和目標檢測算法能監測松材線蟲病的發生現狀和染病松樹的分布情況，為后續的實時處理提供了依據，具有重要的現實意義。

松樹感染松材線蟲病是一個動態過程，不同的染病階段癥狀不盡相同[37-38]，染病初期只有少量松樹針葉褪色發黃，染病中期大部分針葉變為黃褐色，染病后期全部針葉變為黃褐色或紅褐色[18]。本方法根據顏色特征只能識別出中后期的染病松樹，無法檢測顏色變化較小的染病初期的松樹，因此不能全面對染病松樹進行檢測。在今后的研究中將采用高光譜數據進行實驗。

5 結束語

使用無人機遙感技術獲取超高空間分辨率的松林影像，并結合深度學習的目標檢測技術對受害木進行檢測，能有效提高受害木的識別效率，且具有較高的檢測精度。當根據受害木冠幅的實際尺寸修改RPN網絡中的錨框尺寸后，在數據集中加入其他枯死樹和紅色闊葉樹的樣本信息，能提高受害木的識別效果和檢測精度，受害木識別的總體精度達到82.42%。使用坐標轉換的方法對預測結果中的受害木進行定位，結合點位合并過程，最終正確定位出494棵受害木。本文方法能及時發現染病松樹，并確定其分布情況，有效監測松材線蟲病疫情的發展動態，可為松林管理人員和森林防護人員及時提供準確的信息，同時也為松材線蟲病災害損失評估和松林管理部門制定松材線蟲病防控目標等提供客觀依據。