針對糖尿病視網膜病變（DR）分級任務中不同種類之間差異性微小特點，提出一種基于跨層雙線性池化（CHBP）的視網膜病變分級算法。首先根據霍夫圓變換（HCT）對輸入圖像進行裁剪，再使用預處理方法提升圖像對比度；然后以擠壓激勵分組殘差網絡（SEResNeXt）作為模型的主干，引入跨層雙線性池化模塊進行分類；最后在訓練過程中引入隨機拼圖生成器進行漸進訓練，并采用中心損失（CL）和焦點損失（FL）方法進一步提升最終分類效果。實驗結果顯示，本文方法在印度糖尿病視網膜病變圖像數據集（IDRiD）中二次加權卡帕系數（QWK）為90.84%，在梅西多數據集（Messidor-2）中受試者工作特征曲線下的面積（AUC）為88.54%。實驗證明，本文提出的算法在糖尿病視網膜病變分級領域具有一定應用價值。
引用本文： 梁禮明, 彭仁杰, 馮駿, 尹江. 基于跨層雙線性池化的糖尿病視網膜病變分級算法研究. 生物醫學工程學雜志, 2022, 39(5): 928-936. doi: 10.7507/1001-5515.202104038 復制
引言
糖尿病視網膜病變（diabetic retinopathy，DR）是由糖尿病引起視網膜微血管病變的眼部疾病，近期研究表明它可能是導致多數糖尿病患者失明的主要原因[1]，而微動脈瘤、硬滲出物等疾病都與DR密切相關[2]。根據近年相關治療結果分析，對DR開展早期治療能夠有效地預防視力障礙并大幅度減小失明的風險。在臨床實踐中，DR診斷需要眼科醫生對患者進行詳細的眼底檢查，其診斷過程是一項非常耗時的任務。因此利用計算機輔助診斷DR從而提高診斷效率，這一需求顯得越來越迫切[3]。
隨著深度學習技術的發展，已證實卷積神經網絡（convolutional neural network，CNN）能夠有效地學習不同DR分級特征，是計算機輔助診斷的有力工具之一。相對于傳統圖像識別領域以手工設計分類特征(如顏色、形狀、黃斑和血管等)的分級方法，基于CNN的方法能夠顯著地提升DR分類的性能[4-5]。例如，Zhou等[6]采用分類和回歸損失的方法，提出一種多網格多任務視網膜分級框架。該框架采用多個網格結構逐步提高圖像分辨率和網絡的深度，同時引入均方差損失函數最小化標簽值和預測值的差異，在加速訓練的同時提高預測精度。而后Ren等[7]提出一種以矢量量化半監督學習方法用于糖尿病黃斑水腫分級，對疑似眼底滲出區域進行分割并矢量化，再通過半監督學習及圖形分類器復判，最終對滲出物的位置及黃斑判斷病變分級嚴重程度。文獻[8-9]的研究中，先利用CNN進行血管分割，再從血管結構中提取鑒別信息，最后根據鑒別信息對糖尿病視網膜進行病變分級。以上這些針對視網膜病變分級的方法雖采用不同深度的網絡結構及分級策略，但要實現準確診斷的目的，依然存在以下幾點局限：① 與粗粒度分類（coarse-grained classification，CGC）不同，DR分級是細粒度分類（fine-grained classification，FGC）任務，具有較大的類內方差；② 醫學領域帶圖像標注的DR分級開源數據集稀少且標注信息有限；③ 醫學領域公開數據集內的數據種類分布極度不均衡。
針對上述問題，本文提出一種基于跨層雙線性池化（cross-layer hierarchical bilinear pooling，CHBP）的DR-FGC算法。由于臨床中DR嚴重程度的分級依據主要取決于小病變區域，如微動脈瘤、出血等，因此采用CHBP模塊互補輸出，能最大限度減少特征丟失并突出小區域細節特征；同時引入漸進訓練策略，配合隨機拼圖技術逐階段優化網絡各個部分參數，引導網絡從細粒度逐步到粗粒度學習特征。為了進一步提升DR分級效果，本文還融合中心損失（center loss，CL）函數和焦點損失（focal loss，FL）函數分別緩解數據集類內方差大和數據類別分布極度不均衡等問題。相對其他現有算法，通過在兩個公開數據集中對本文算法進行驗證，以期實現在臨床診斷中對DR快速篩查的目的。
1 數據來源及預處理
1.1 數據來源
本文算法使用2018年國際生物醫學成像技術研討會（international symposium on biomedical imaging，ISBI）上印度DR圖像挑戰賽開源數據集（Indian DR image dataset，IDRiD）（網址為：
圖1
不同嚴重等級的DR圖
Figure1.
Diagrams of different severity levels of DR
IDRiD數據集中訓練集圖像共413張，其中健康圖像134張、輕度非增殖性圖像20張、中度非增殖性圖像136張、重度非增殖性圖像74張、增殖性圖像49張。測試集圖像共103張，其中健康圖像34張、輕度非增殖性圖像5張、中度非增殖性圖像32張、重度非增殖性圖像19張、增殖性圖像13張。Messidor-2數據集由來自于三個眼科部門的1 200張眼底彩色圖像組成，其分辨率分別為1 440 × 960，2 240 × 1 488和2 304 × 1 536。
DR根據微動脈瘤數量（μA）、出血次數（H）、新生血管（NV=1）和無新生血管（NV=0）情況，將上述圖像分為4級，其中第0級健康圖像（μA=0和H=0）共546張、第1級圖像（0<μA≤5和H=0）共254張、第2級圖像（5<μA≤15或0<H<5和NV=0）共247張、第3級圖像（μA≥15或H≥5和NV=1）共153張，對比IDRiD數據集并未出現嚴重類別不均衡現象。在DR篩查過程中，第0級健康圖像和第1級圖像是計算機輔助設計系統和臨床專家最難區分的任務。Sánchez等[13]將第0級健康圖像和第1級圖像作為可參考圖像分為一組，把第2級和第3級圖像作為不可參考分為另一組。這兩類設定已經廣泛應用在現有DR篩選中，因此本文將Messidor-2數據集作為二分類任務[14]。
1.2 數據預處理算法原理
針對DR數據集中存在大量黑色背景問題，首先使用霍夫圓變換（hough circle transform，HCT）定位圖像中眼球位置并裁剪多余的黑色背景區域及部分眼球邊緣區域像素，再將圖像分辨率全部調整到448 × 448。同時為了突出視網膜中小病變區域特征，對圖像I進行高斯濾波操作以平滑中小病變區域，得到的高斯模糊圖作為眼球平均背景圖像，然后將加權的原始圖像減去加權高斯模糊圖即可獲得眼球特征圖，最后為了避免產生負值引入常數γ，如式(1)～式(2)所示：
 |
 |
式中，α、β和γ分別代表線性加權值，其值分別取4、?4和128；*表示濾波操作；x、y分別代表圖像單個像素點x軸和y軸坐標，σ為標準差值，e為自然常數，Gσ為二維高斯核，
為經過加權融合后的圖像，如圖2所示。
圖2
圖像預處理
Figure2.
Image preprocessing
2 網絡體系結構
2.1 主干網絡原理
對比其他圖像分類任務，DR分級由于各類圖像特征差異小、數據集樣本量少和類別分布不均勻等問題變得更具挑戰性。相比以殘差網絡（residual network，ResNet）作為主干網絡，擠壓激勵分組ResNet（squeeze-and-excitation grouping ResNet，SEResNeXt）不但融合了通道注意力機制還增大了網絡的基數，使模型辨別能力到達一個更高的水平[15-16]。在SEResNeXt瓶頸結構中，繼承ResNet經典的堆砌結構同時，運用分組卷積將原有的3 × 3卷積操作拆分成數個部分。具體而言，上述操作分為拆分、轉換和聚合三部分：首先設定需要拆分的向量個數D，將向量x拆分成D個低維向量 
，然后將多個低維度向量經不同的卷積 
提取特征得到 
，最后把多個低維度向量特征按通道方向連接，并由1 × 1卷積 
將各個低維向量融合形成 
，如式(3)所示：
 |
相比于ResNet的瓶頸結構，SEResNeXt在大幅度精簡網絡參數量的同時顯著提升網絡識別的準確率。而擠壓激勵通道注意力模塊（squeeze-and-excitation channel attention module，SECAM）通過提取特征通道之間的相互依賴關系，允許網絡對通道特征進行重新校準，通過全局信息來加強有用的特征并抑制不太重要的特征。給定特征輸入矩陣為u，將其進行全局平均池化操作獲取全局特征信息向量，如式(4)所示：
 |
式中，i和j分別代表特征輸入u中對應坐標點的像素值，H和W代表特征輸入矩陣的長和寬，c代表特征輸入對應的通道維度，最終特征輸出為包含通道信息的 
。為了提取通道中的依賴關系，使用兩層全連接層對通道特征 
進行學習，最后用S型生長曲線(sigmoid)激活函數作為門控機制，如式(5)所示：
 |
式中，δ代表線性整流激活函數（linear rectification function，ReLU），
和 
分別為全連接層學習參數矩陣，同時為了限制模型復雜度將輸入的向量維度c減少至r倍，并由第二個全連接層恢復至c個通道維度，σ代表sigmoid激活函數，最終得到歸一化后的通道注意力特征矩陣s。將通道注意力特征矩陣s與特征u逐通道相乘，獲得最終調整通道后的特征圖。網絡瓶頸層結構圖如圖3所示。
圖3
SEResNeXt主干網絡瓶頸層結構圖
Figure3.
Structure of SEResNeXt backbone network bottleneck layer
2.2 跨層雙線性池化模塊
以往的分類CNN通常在主干網絡后面再添加全局平均池化層降維并通過全連接層直接對特征進行分類，而作為小樣本細粒度的DR分級任務顯然不適用。首先，僅將最后一層卷積的激活值作為圖像的特征表示并采用池化層對特征進行縮減，會導致特征大量丟失，不足以描述對象細粒度的各種語義表征；其次，它忽略了層間部分特征交互關系和細粒度特征學習的相互關聯性。中間層卷積激活的語義信息可以通過互補的方式形成更具判別的特征，這部分信息對于細粒度視覺識別具有重要意義。因此本文采用一種CHBP方法來捕獲層間部分特征關系，這與其他分類方法相比具有更優的性能[17]。
為了最大限度地減少細粒度識別信息的損失，CHBP框架集成了多個跨層雙線性特征，以增強它們的表示能力。由于來自不同卷積層的特征具有互補性，充分利用中間卷積層的激活特征值有助于鑒別特征學習，因此所有交互層的雙線性特征在最終分類之前被級聯。三條不同的特征通過獨立的線性映射，并在層間相互交互，強化相同位置的共有特征并抑制無關特征，形成自我注意力機制。本文所提出的網絡受益于層間特征交互和細粒度特征學習之間的相互增強。
跨層雙線性模塊融合SEResNext網絡最后三個瓶頸層的激活函數特征，其中雙線性向量的輸出連接如式(6)所示：
 |
式中，
、
、
分別代表輸入的特征圖；V、U、S是將特征映射成雙線性向量的投影矩陣；P是分類參數矩陣；
為哈達瑪積；concat(·)代表特征圖按照通道維度連接操作；Z為跨層雙線性模塊最終輸出結果。
2.3 損失函數
在深度學習損失函數選擇方面，由于DR數據集自身存在嚴重的數據類別不均衡問題，同時每類DR圖像之間特征差異性小且類別粒度相對于其他分類任務更為精細，使得傳統交叉熵損失函數訓練的網絡難以區分DR每一種類別。為了解決上述問題，本文提出一種新的組合損失函數，即CL函數和FL函數。CL函數是為了提高CNN判別特征空間的鑒別能力而提出的[18]。具體來說，CL函數能夠同時學習每個類深度特征中心，并懲罰深度特征中心與其對應的類中心的距離。從而迫使不同類的深度特征保持分離，進而擴大類間特征的差異性，并將類內樣本之間的距離拉近，使得網絡學習特征的鑒別能力能夠高度增強，如式(7)所示：
 |
式中，m對應分類任務類別總數，
代表種類為i的深層特征的類別中心，用于表征類內變化，
代表CNN預測對應種類i的類別中心，
為CL值。
在網絡優化過程中，
隨著深度特征的變化而更新，但由于類中心 
代表所有樣本特征的平均值，理論上無法直接獲取，因此在每個批次中計算當前數據與中心的距離，然后以梯度的形式加到中心上，類似于梯度下降的參數修正，同時再增加一個標量α用于避免小批量數據帶來劇烈抖動，x代表網絡總體預測值，
為 
求偏導，
為x求偏導，如式(8)所示：
 |
由于數據集中存在嚴重數據類別不均衡問題，FL損失函數通過加權誤差調制指數來重塑交叉熵，以減輕不均衡所帶來的干擾[19]，如式(9)所示：
 |
式中，γ代表加權誤差調制系數，
代表對應類別的標簽，
代表CNN的預測值，
為FL損失函數結果。
最終損失函數是度量CL函數和FL函數的加權組合，引入超參數λ以均衡兩個損失函數的值，如式(10)所示：
 |
式中，Ltotal為最終的融合損失，超參數設置值λ為0.001。
2.4 漸進訓練策略
由于DR分級是一個漸進的過程，不同種類之間存在微妙的變化關系，而同屬一類的眼球圖像在形態也有很大差異。作為FGC任務，近年研究主要集中在定位更具判別性、互補和各種尺寸的特征上。相比于首先定位差異性特征的方式，本文算法側重發掘圖像中分布的大尺度結構，并將注意轉移到越來越小的尺寸信息中，引導網絡從細粒度逐步到粗粒度學習特征。
漸進訓練策略最早運用于生成對抗網絡，從低階段網絡開始訓練，然后通過添加卷積層逐步訓練更高階段網絡，相比以往的訓練方式能夠通過中間監督簡化網絡中信息傳遞[20]。由于低階段網絡的感受野和表達能力有限，網絡被迫從局部細節中尋找可利用的差異性信息，這種漸進訓練策略允許模型將鑒別信息從局部細節定位到全局結構中。漸進訓練策略將整個訓練過程分為4個階段，具體流程如圖4所示。階段1首先訓練步驟1中的低層次卷積層，其中兩層卷積堆疊的卷積塊（convolution block，CB）用于通道特征融合和降維，然后使用全局最大池化層和全連接層分類器對步驟1的網絡參數訓練更新；然后按照上述操作依次訓練步驟2和步驟3中的卷積層并同步更新權重，損失函數1、損失函數2和損失函數3均為交叉熵損失；最后訓練階段4中步驟4所對應的卷積層參數，包含整個主干網絡和雙線性池化分類器，最終損失包含FL和CL的組合損失函數。在整個訓練過程中為了進一步提升效果，引入隨機拼圖生成模塊，迫使模型在每個訓練階段學習特定粒度的信息。
圖4
漸進訓練策略流程
Figure4.
Progressive training strategy process
2.5 隨機拼圖生成模塊
圖像拼圖處理最早運用于不同漸進訓練方式的自監督任務[21]。Chen等[22]在精細圖像分類任務中采用破壞全局結構的圖像切片操作，雖取得了較好的分類效果，但在訓練過程中切片操作意味著很難獲取多種尺度的區域特征。為了避免上述情況，本文在漸進訓練過程中引入隨機拼圖生成器，用于限制網絡學習特定區域信息并加強局部細節特征和圖像重建，以學習局部區域之間的語義相關性。首先，隨機拼圖生成器將輸入圖像分割成數個小塊用于訓練模型底層網絡，然后逐步減小拼圖分割數量并開始對應高層網絡訓練，其目的是通過不同粒度區域迫使模型在每個訓練步驟中學習特定粒度的信息。具體操作如下：給定一張輸入圖像，將其按照長寬方向平均分割，形成多個相同大小的小塊，然后將每個小塊的位置進行隨機洗牌合并成一個全新圖像。關于每個階段小塊數量選擇遵守下面兩個條件：① 切塊的大小應小于相應階段的感受野，否則會降低拼圖生成器性能。② 切片大小應隨各階段感受野的增加而同比例增加。通常每個階段感受野大概是上一階段的兩倍，故設定數量分別為16、4、1。
在推理階段網絡中不再使用拼圖生成和漸進學習策略，直接通過主干網絡和雙線性分類器輸出最終結果。
2.6 跨層雙線性池化的糖尿病視網膜病變分級網絡
本文提出的CHBP網絡由兩部分組成：SEResNeXt主干網絡和CHBP分類器。其中SEResNeXt主干網絡采用免費公開的圖像網絡數據集（ImageNet dataset）（網址為：
圖5
網絡大體結構
Figure5.
General structure of the network
3 實驗內容與結果分析
3.1 評價指標
為了準確地評估本文所提方法在IDRiD和Messidor-2數據集上的表現，同時方便與之前的方法對比，本文對上述兩個數據集進行5折交叉驗證。針對不同數據集評判指標的差異性，IDRiD數據集利用二次加權卡帕系數（quadratic weighted Kappa，QWK）進行一致性檢驗評估，其范圍在–1～1之間。通過其衡量兩個評級之間的一致性，值越大，一致性越高，如式(11)～式(12)所示：
 |
 |
式中，N代表分類的總類數，i和j代表預測類別和標簽類別，
代表每一對(i, j)的懲罰權重，
代表第i類判別為j類的個數，
代表i類總數乘以j類總數再除以總個數。
同時在Messidor-2數據集上評估時采用靈敏度（sensitivity，Sen）、特異性（specificity，Spe）、真陽性率（true positive rate，TPR）、假陽性率（false positive rate，FPR）和受試者工作特征曲線（receiver operating characteristic curve，ROC）下的面積（area under curve，AUC）值作為指標評價分類性能，對模型進行評判。如式(13)～式(15)所示：
 |
 |
 |
式中，真陽性（true positive，TP）代表正樣本被正確分類的個數，真陰性（true negative，TN）代表負樣本被正確分類的個數，假陽性（false positive，FP）代表正樣本被錯誤分類的個數，假陰性（false negative，FN）代表負樣本被錯誤分類的個數。通過以TPR為縱軸，以FPR為橫軸繪制得到ROC曲線，求出ROC曲線下面積獲得AUC值。
3.2 訓練細節
針對數據集中圖像分辨率大小不一致的問題，為了方便訓練和測試，本文將其統一調整分辨率為448×448。同時在訓練過程中使用數據增強的手段，包括：隨機調整大小裁剪、隨機仿射、隨機水平翻轉、隨機垂直翻轉和高斯模糊等。本文所有的實驗均在操作系統Ubuntu 16.04（Canonical Inc.，美國）進行；建模基于深度學習架構Pytorch 1.5（Facebook Inc.，美國）和計算統一設備架構CUDA 10.1（Nvidia Inc.，美國）。計算機具體配置：顯卡（Nvidia GeFore GTX2070 GPU，Nvidia Inc.，美國）、中央處理器（Intel Core TM i7-6700H CPU, Inter Inc.，美國）。所有骨干網絡均在ImageNet dataset數據集上進行預先訓練。對于優化器采用自適應矩估計（adaptive moment estimation，Adam），動量（momentum，Mom）設置為0.9。批量大小和輪次分別設置為4和200。學習率初始值為0.01，并采用循環學習率訓練策略[23]。
3.3 實驗結果分析
為了更直觀地分析CHBP算法中不同模塊對分級性能的貢獻，本文設計了四組對照實驗：① 僅僅保留CHBP算法主干網絡（實驗1）；② 將CHBP算法的組合損失函數改成交叉熵損失函數，并去掉漸進訓練方法（實驗2）；③ 只去掉CHBP算法的漸進訓練方法（實驗3）；④ 本文所提完整的CHBP算法（實驗4）。分別對四種不同的對照實驗進行仿真，結果如表1所示。
表1
對照實驗在IDRiD數據集中5折交叉驗證表現
Table1.
Five-fold cross-validation performance of control experiment in the IDRiD dataset
表1反映了CHBP算法中不同模塊在IDRiD數據集上性能貢獻情況，其中由實驗1和實驗2可知雙線性池化層能夠使QWK提升1.89%，在所有模塊中提升最明顯。而根據實驗2和實驗3可以看出FL函數和CL函數能夠有效地提升DR性能。最后引入漸進訓練方式引導網絡學習，對于最終的結果也有一定的貢獻。
3.4 對比分析
為了充分說明本文算法對于DR分級性能，表2和表3分別給出本文模型與其他模型在IDRiD和Messidor-2數據集上實驗結果對比。
表2
不同模型在IDRiD數據集中5折交叉驗證表現
Table2.
Five-fold cross-validation performance of different models in the IDRiD dataset
表3
不同模型在Messidor-2數據集中的表現
Table3.
Performance of different models in the Messidor-2 dataset
如表2所示，本文在IDRiD數據集上對比現有主流FGC算法和最新DR分級算法的分類結果，其QWK均低于本文所提出的算法。文獻[17]提出雙線性池化模塊對網絡輸出特征進行融合，但由于自身主干網絡特征提取能力有限，而與本文算法性能差異最大。文獻[24]提出一種利用局部信息和全局外觀信息的非對稱多分支結構并引入弱監督增強學習能力，但是無法將類間特征分離且無弱監督能力，雖然有較大進步但實驗結果相比于本文算法下降了0.028 0。文獻[25]提出漸進訓練方法引導網絡學習，但特征融合時層間特征差距過大且部分特征未經充分提取也使得最終精度受到影響。文獻[26]采用多主干模型增強分類效果，最終結果依然弱于本文所提算法。文獻[27]融合多種FGC策略并通過度量學習獲得更具鑒別力的特征以幫助最終分類。對比上述研究，由分類結果可知，本文方法在IDRiD數據集中取得更有競爭力的性能結果。5折交叉驗證（交叉驗證-1～交叉驗證-5）訓練圖及學習率變化趨勢，如圖6所示。
圖6
學習率變化圖及模型交叉驗證QWK值變化趨勢圖
Figure6.
Learning rate change and model cross-validation QWK value change trend
如表3所示，現有FGC算法和最新的DR分級算法在Messidor-2數據集中，本文Sen對比文獻[28]和文獻[29]分別提高0.108 0和0.202 3。但是特異性低于文獻[29]，AUC均高于其他兩種算法，分別增長0.032 4和0.006 7。這也表明本文算法在Messidor-2數據集上依然具有一定的泛化能力。
4 結束語
本文提出了一個新的CHBP的DR-FGC模型，其整體結構由主干網絡（SEResNeXt）和CHBP分類器組成，并引入CL函數和FL函數緩解數據不均衡和種類差異度小的問題，同時利用漸進訓練方法提升精度。實驗結果證明，CHBP算法能夠有效地根據DR不同種類之間微小差異判斷DR等級。由于本文網絡整體參數量大，部分網絡還有精減的空間，而且引入漸進訓練方法大幅度增加了訓練時間，因此下一步研究重心將集中在精簡網絡和加快推理上。本文模型不但提高了DR分級的效率，而且不會出現醫生由于主觀原因而導致誤判，在DR分級領域具有很大的前景。
重要聲明
利益沖突聲明：本文全體作者均聲明不存在利益沖突。
作者貢獻聲明：梁禮明主要負責項目主持、平臺搭建、算法程序設計、協調溝通以及計劃安排，彭仁杰主要負責實驗流程、數據記錄與分析、論文編寫以及算法程序設計；馮駿和尹江主要負責提供實驗指導，數據分析指導，論文審閱修訂。
引言
糖尿病視網膜病變（diabetic retinopathy，DR）是由糖尿病引起視網膜微血管病變的眼部疾病，近期研究表明它可能是導致多數糖尿病患者失明的主要原因[1]，而微動脈瘤、硬滲出物等疾病都與DR密切相關[2]。根據近年相關治療結果分析，對DR開展早期治療能夠有效地預防視力障礙并大幅度減小失明的風險。在臨床實踐中，DR診斷需要眼科醫生對患者進行詳細的眼底檢查，其診斷過程是一項非常耗時的任務。因此利用計算機輔助診斷DR從而提高診斷效率，這一需求顯得越來越迫切[3]。
隨著深度學習技術的發展，已證實卷積神經網絡（convolutional neural network，CNN）能夠有效地學習不同DR分級特征，是計算機輔助診斷的有力工具之一。相對于傳統圖像識別領域以手工設計分類特征(如顏色、形狀、黃斑和血管等)的分級方法，基于CNN的方法能夠顯著地提升DR分類的性能[4-5]。例如，Zhou等[6]采用分類和回歸損失的方法，提出一種多網格多任務視網膜分級框架。該框架采用多個網格結構逐步提高圖像分辨率和網絡的深度，同時引入均方差損失函數最小化標簽值和預測值的差異，在加速訓練的同時提高預測精度。而后Ren等[7]提出一種以矢量量化半監督學習方法用于糖尿病黃斑水腫分級，對疑似眼底滲出區域進行分割并矢量化，再通過半監督學習及圖形分類器復判，最終對滲出物的位置及黃斑判斷病變分級嚴重程度。文獻[8-9]的研究中，先利用CNN進行血管分割，再從血管結構中提取鑒別信息，最后根據鑒別信息對糖尿病視網膜進行病變分級。以上這些針對視網膜病變分級的方法雖采用不同深度的網絡結構及分級策略，但要實現準確診斷的目的，依然存在以下幾點局限：① 與粗粒度分類（coarse-grained classification，CGC）不同，DR分級是細粒度分類（fine-grained classification，FGC）任務，具有較大的類內方差；② 醫學領域帶圖像標注的DR分級開源數據集稀少且標注信息有限；③ 醫學領域公開數據集內的數據種類分布極度不均衡。
針對上述問題，本文提出一種基于跨層雙線性池化（cross-layer hierarchical bilinear pooling，CHBP）的DR-FGC算法。由于臨床中DR嚴重程度的分級依據主要取決于小病變區域，如微動脈瘤、出血等，因此采用CHBP模塊互補輸出，能最大限度減少特征丟失并突出小區域細節特征；同時引入漸進訓練策略，配合隨機拼圖技術逐階段優化網絡各個部分參數，引導網絡從細粒度逐步到粗粒度學習特征。為了進一步提升DR分級效果，本文還融合中心損失（center loss，CL）函數和焦點損失（focal loss，FL）函數分別緩解數據集類內方差大和數據類別分布極度不均衡等問題。相對其他現有算法，通過在兩個公開數據集中對本文算法進行驗證，以期實現在臨床診斷中對DR快速篩查的目的。
1 數據來源及預處理
1.1 數據來源
本文算法使用2018年國際生物醫學成像技術研討會（international symposium on biomedical imaging，ISBI）上印度DR圖像挑戰賽開源數據集（Indian DR image dataset，IDRiD）（網址為：
圖1
不同嚴重等級的DR圖
Figure1.
Diagrams of different severity levels of DR
IDRiD數據集中訓練集圖像共413張，其中健康圖像134張、輕度非增殖性圖像20張、中度非增殖性圖像136張、重度非增殖性圖像74張、增殖性圖像49張。測試集圖像共103張，其中健康圖像34張、輕度非增殖性圖像5張、中度非增殖性圖像32張、重度非增殖性圖像19張、增殖性圖像13張。Messidor-2數據集由來自于三個眼科部門的1 200張眼底彩色圖像組成，其分辨率分別為1 440 × 960，2 240 × 1 488和2 304 × 1 536。
DR根據微動脈瘤數量（μA）、出血次數（H）、新生血管（NV=1）和無新生血管（NV=0）情況，將上述圖像分為4級，其中第0級健康圖像（μA=0和H=0）共546張、第1級圖像（0<μA≤5和H=0）共254張、第2級圖像（5<μA≤15或0<H<5和NV=0）共247張、第3級圖像（μA≥15或H≥5和NV=1）共153張，對比IDRiD數據集并未出現嚴重類別不均衡現象。在DR篩查過程中，第0級健康圖像和第1級圖像是計算機輔助設計系統和臨床專家最難區分的任務。Sánchez等[13]將第0級健康圖像和第1級圖像作為可參考圖像分為一組，把第2級和第3級圖像作為不可參考分為另一組。這兩類設定已經廣泛應用在現有DR篩選中，因此本文將Messidor-2數據集作為二分類任務[14]。
1.2 數據預處理算法原理
針對DR數據集中存在大量黑色背景問題，首先使用霍夫圓變換（hough circle transform，HCT）定位圖像中眼球位置并裁剪多余的黑色背景區域及部分眼球邊緣區域像素，再將圖像分辨率全部調整到448 × 448。同時為了突出視網膜中小病變區域特征，對圖像I進行高斯濾波操作以平滑中小病變區域，得到的高斯模糊圖作為眼球平均背景圖像，然后將加權的原始圖像減去加權高斯模糊圖即可獲得眼球特征圖，最后為了避免產生負值引入常數γ，如式(1)～式(2)所示：
|
|
式中，α、β和γ分別代表線性加權值，其值分別取4、?4和128；*表示濾波操作；x、y分別代表圖像單個像素點x軸和y軸坐標，σ為標準差值，e為自然常數，Gσ為二維高斯核，為經過加權融合后的圖像，如圖2所示。
圖2
圖像預處理
Figure2.
Image preprocessing
2 網絡體系結構
2.1 主干網絡原理
對比其他圖像分類任務，DR分級由于各類圖像特征差異小、數據集樣本量少和類別分布不均勻等問題變得更具挑戰性。相比以殘差網絡（residual network，ResNet）作為主干網絡，擠壓激勵分組ResNet（squeeze-and-excitation grouping ResNet，SEResNeXt）不但融合了通道注意力機制還增大了網絡的基數，使模型辨別能力到達一個更高的水平[15-16]。在SEResNeXt瓶頸結構中，繼承ResNet經典的堆砌結構同時，運用分組卷積將原有的3 × 3卷積操作拆分成數個部分。具體而言，上述操作分為拆分、轉換和聚合三部分：首先設定需要拆分的向量個數D，將向量x拆分成D個低維向量 ，然后將多個低維度向量經不同的卷積 提取特征得到 ，最后把多個低維度向量特征按通道方向連接，并由1 × 1卷積 將各個低維向量融合形成 ，如式(3)所示：
|
相比于ResNet的瓶頸結構，SEResNeXt在大幅度精簡網絡參數量的同時顯著提升網絡識別的準確率。而擠壓激勵通道注意力模塊（squeeze-and-excitation channel attention module，SECAM）通過提取特征通道之間的相互依賴關系，允許網絡對通道特征進行重新校準，通過全局信息來加強有用的特征并抑制不太重要的特征。給定特征輸入矩陣為u，將其進行全局平均池化操作獲取全局特征信息向量，如式(4)所示：
|
式中，i和j分別代表特征輸入u中對應坐標點的像素值，H和W代表特征輸入矩陣的長和寬，c代表特征輸入對應的通道維度，最終特征輸出為包含通道信息的 。為了提取通道中的依賴關系，使用兩層全連接層對通道特征 進行學習，最后用S型生長曲線(sigmoid)激活函數作為門控機制，如式(5)所示：
|
式中，δ代表線性整流激活函數（linear rectification function，ReLU）， 和 分別為全連接層學習參數矩陣，同時為了限制模型復雜度將輸入的向量維度c減少至r倍，并由第二個全連接層恢復至c個通道維度，σ代表sigmoid激活函數，最終得到歸一化后的通道注意力特征矩陣s。將通道注意力特征矩陣s與特征u逐通道相乘，獲得最終調整通道后的特征圖。網絡瓶頸層結構圖如圖3所示。
圖3
SEResNeXt主干網絡瓶頸層結構圖
Figure3.
Structure of SEResNeXt backbone network bottleneck layer
2.2 跨層雙線性池化模塊
以往的分類CNN通常在主干網絡后面再添加全局平均池化層降維并通過全連接層直接對特征進行分類，而作為小樣本細粒度的DR分級任務顯然不適用。首先，僅將最后一層卷積的激活值作為圖像的特征表示并采用池化層對特征進行縮減，會導致特征大量丟失，不足以描述對象細粒度的各種語義表征；其次，它忽略了層間部分特征交互關系和細粒度特征學習的相互關聯性。中間層卷積激活的語義信息可以通過互補的方式形成更具判別的特征，這部分信息對于細粒度視覺識別具有重要意義。因此本文采用一種CHBP方法來捕獲層間部分特征關系，這與其他分類方法相比具有更優的性能[17]。
為了最大限度地減少細粒度識別信息的損失，CHBP框架集成了多個跨層雙線性特征，以增強它們的表示能力。由于來自不同卷積層的特征具有互補性，充分利用中間卷積層的激活特征值有助于鑒別特征學習，因此所有交互層的雙線性特征在最終分類之前被級聯。三條不同的特征通過獨立的線性映射，并在層間相互交互，強化相同位置的共有特征并抑制無關特征，形成自我注意力機制。本文所提出的網絡受益于層間特征交互和細粒度特征學習之間的相互增強。
跨層雙線性模塊融合SEResNext網絡最后三個瓶頸層的激活函數特征，其中雙線性向量的輸出連接如式(6)所示：
|
式中，、、 分別代表輸入的特征圖；V、U、S是將特征映射成雙線性向量的投影矩陣；P是分類參數矩陣；為哈達瑪積；concat(·)代表特征圖按照通道維度連接操作；Z為跨層雙線性模塊最終輸出結果。
2.3 損失函數
在深度學習損失函數選擇方面，由于DR數據集自身存在嚴重的數據類別不均衡問題，同時每類DR圖像之間特征差異性小且類別粒度相對于其他分類任務更為精細，使得傳統交叉熵損失函數訓練的網絡難以區分DR每一種類別。為了解決上述問題，本文提出一種新的組合損失函數，即CL函數和FL函數。CL函數是為了提高CNN判別特征空間的鑒別能力而提出的[18]。具體來說，CL函數能夠同時學習每個類深度特征中心，并懲罰深度特征中心與其對應的類中心的距離。從而迫使不同類的深度特征保持分離，進而擴大類間特征的差異性，并將類內樣本之間的距離拉近，使得網絡學習特征的鑒別能力能夠高度增強，如式(7)所示：
|
式中，m對應分類任務類別總數， 代表種類為i的深層特征的類別中心，用于表征類內變化， 代表CNN預測對應種類i的類別中心， 為CL值。 在網絡優化過程中， 隨著深度特征的變化而更新，但由于類中心 代表所有樣本特征的平均值，理論上無法直接獲取，因此在每個批次中計算當前數據與中心的距離，然后以梯度的形式加到中心上，類似于梯度下降的參數修正，同時再增加一個標量α用于避免小批量數據帶來劇烈抖動，x代表網絡總體預測值， 為 求偏導， 為x求偏導，如式(8)所示：
|
由于數據集中存在嚴重數據類別不均衡問題，FL損失函數通過加權誤差調制指數來重塑交叉熵，以減輕不均衡所帶來的干擾[19]，如式(9)所示：
|
式中，γ代表加權誤差調制系數， 代表對應類別的標簽， 代表CNN的預測值， 為FL損失函數結果。
最終損失函數是度量CL函數和FL函數的加權組合，引入超參數λ以均衡兩個損失函數的值，如式(10)所示：
|
式中，Ltotal為最終的融合損失，超參數設置值λ為0.001。
2.4 漸進訓練策略
由于DR分級是一個漸進的過程，不同種類之間存在微妙的變化關系，而同屬一類的眼球圖像在形態也有很大差異。作為FGC任務，近年研究主要集中在定位更具判別性、互補和各種尺寸的特征上。相比于首先定位差異性特征的方式，本文算法側重發掘圖像中分布的大尺度結構，并將注意轉移到越來越小的尺寸信息中，引導網絡從細粒度逐步到粗粒度學習特征。
漸進訓練策略最早運用于生成對抗網絡，從低階段網絡開始訓練，然后通過添加卷積層逐步訓練更高階段網絡，相比以往的訓練方式能夠通過中間監督簡化網絡中信息傳遞[20]。由于低階段網絡的感受野和表達能力有限，網絡被迫從局部細節中尋找可利用的差異性信息，這種漸進訓練策略允許模型將鑒別信息從局部細節定位到全局結構中。漸進訓練策略將整個訓練過程分為4個階段，具體流程如圖4所示。階段1首先訓練步驟1中的低層次卷積層，其中兩層卷積堆疊的卷積塊（convolution block，CB）用于通道特征融合和降維，然后使用全局最大池化層和全連接層分類器對步驟1的網絡參數訓練更新；然后按照上述操作依次訓練步驟2和步驟3中的卷積層并同步更新權重，損失函數1、損失函數2和損失函數3均為交叉熵損失；最后訓練階段4中步驟4所對應的卷積層參數，包含整個主干網絡和雙線性池化分類器，最終損失包含FL和CL的組合損失函數。在整個訓練過程中為了進一步提升效果，引入隨機拼圖生成模塊，迫使模型在每個訓練階段學習特定粒度的信息。
圖4
漸進訓練策略流程
Figure4.
Progressive training strategy process
2.5 隨機拼圖生成模塊
圖像拼圖處理最早運用于不同漸進訓練方式的自監督任務[21]。Chen等[22]在精細圖像分類任務中采用破壞全局結構的圖像切片操作，雖取得了較好的分類效果，但在訓練過程中切片操作意味著很難獲取多種尺度的區域特征。為了避免上述情況，本文在漸進訓練過程中引入隨機拼圖生成器，用于限制網絡學習特定區域信息并加強局部細節特征和圖像重建，以學習局部區域之間的語義相關性。首先，隨機拼圖生成器將輸入圖像分割成數個小塊用于訓練模型底層網絡，然后逐步減小拼圖分割數量并開始對應高層網絡訓練，其目的是通過不同粒度區域迫使模型在每個訓練步驟中學習特定粒度的信息。具體操作如下：給定一張輸入圖像，將其按照長寬方向平均分割，形成多個相同大小的小塊，然后將每個小塊的位置進行隨機洗牌合并成一個全新圖像。關于每個階段小塊數量選擇遵守下面兩個條件：① 切塊的大小應小于相應階段的感受野，否則會降低拼圖生成器性能。② 切片大小應隨各階段感受野的增加而同比例增加。通常每個階段感受野大概是上一階段的兩倍，故設定數量分別為16、4、1。
在推理階段網絡中不再使用拼圖生成和漸進學習策略，直接通過主干網絡和雙線性分類器輸出最終結果。
2.6 跨層雙線性池化的糖尿病視網膜病變分級網絡
本文提出的CHBP網絡由兩部分組成：SEResNeXt主干網絡和CHBP分類器。其中SEResNeXt主干網絡采用免費公開的圖像網絡數據集（ImageNet dataset）（網址為：
圖5
網絡大體結構
Figure5.
General structure of the network
3 實驗內容與結果分析
3.1 評價指標
為了準確地評估本文所提方法在IDRiD和Messidor-2數據集上的表現，同時方便與之前的方法對比，本文對上述兩個數據集進行5折交叉驗證。針對不同數據集評判指標的差異性，IDRiD數據集利用二次加權卡帕系數（quadratic weighted Kappa，QWK）進行一致性檢驗評估，其范圍在–1～1之間。通過其衡量兩個評級之間的一致性，值越大，一致性越高，如式(11)～式(12)所示：
|
|
式中，N代表分類的總類數，i和j代表預測類別和標簽類別， 代表每一對(i, j)的懲罰權重， 代表第i類判別為j類的個數， 代表i類總數乘以j類總數再除以總個數。
同時在Messidor-2數據集上評估時采用靈敏度（sensitivity，Sen）、特異性（specificity，Spe）、真陽性率（true positive rate，TPR）、假陽性率（false positive rate，FPR）和受試者工作特征曲線（receiver operating characteristic curve，ROC）下的面積（area under curve，AUC）值作為指標評價分類性能，對模型進行評判。如式(13)～式(15)所示：
|
|
|
式中，真陽性（true positive，TP）代表正樣本被正確分類的個數，真陰性（true negative，TN）代表負樣本被正確分類的個數，假陽性（false positive，FP）代表正樣本被錯誤分類的個數，假陰性（false negative，FN）代表負樣本被錯誤分類的個數。通過以TPR為縱軸，以FPR為橫軸繪制得到ROC曲線，求出ROC曲線下面積獲得AUC值。
3.2 訓練細節
針對數據集中圖像分辨率大小不一致的問題，為了方便訓練和測試，本文將其統一調整分辨率為448×448。同時在訓練過程中使用數據增強的手段，包括：隨機調整大小裁剪、隨機仿射、隨機水平翻轉、隨機垂直翻轉和高斯模糊等。本文所有的實驗均在操作系統Ubuntu 16.04（Canonical Inc.，美國）進行；建模基于深度學習架構Pytorch 1.5（Facebook Inc.，美國）和計算統一設備架構CUDA 10.1（Nvidia Inc.，美國）。計算機具體配置：顯卡（Nvidia GeFore GTX2070 GPU，Nvidia Inc.，美國）、中央處理器（Intel Core TM i7-6700H CPU, Inter Inc.，美國）。所有骨干網絡均在ImageNet dataset數據集上進行預先訓練。對于優化器采用自適應矩估計（adaptive moment estimation，Adam），動量（momentum，Mom）設置為0.9。批量大小和輪次分別設置為4和200。學習率初始值為0.01，并采用循環學習率訓練策略[23]。
3.3 實驗結果分析
為了更直觀地分析CHBP算法中不同模塊對分級性能的貢獻，本文設計了四組對照實驗：① 僅僅保留CHBP算法主干網絡（實驗1）；② 將CHBP算法的組合損失函數改成交叉熵損失函數，并去掉漸進訓練方法（實驗2）；③ 只去掉CHBP算法的漸進訓練方法（實驗3）；④ 本文所提完整的CHBP算法（實驗4）。分別對四種不同的對照實驗進行仿真，結果如表1所示。
表1
對照實驗在IDRiD數據集中5折交叉驗證表現
Table1.
Five-fold cross-validation performance of control experiment in the IDRiD dataset
表1反映了CHBP算法中不同模塊在IDRiD數據集上性能貢獻情況，其中由實驗1和實驗2可知雙線性池化層能夠使QWK提升1.89%，在所有模塊中提升最明顯。而根據實驗2和實驗3可以看出FL函數和CL函數能夠有效地提升DR性能。最后引入漸進訓練方式引導網絡學習，對于最終的結果也有一定的貢獻。
3.4 對比分析
為了充分說明本文算法對于DR分級性能，表2和表3分別給出本文模型與其他模型在IDRiD和Messidor-2數據集上實驗結果對比。
表2
不同模型在IDRiD數據集中5折交叉驗證表現
Table2.
Five-fold cross-validation performance of different models in the IDRiD dataset
表3
不同模型在Messidor-2數據集中的表現
Table3.
Performance of different models in the Messidor-2 dataset
如表2所示，本文在IDRiD數據集上對比現有主流FGC算法和最新DR分級算法的分類結果，其QWK均低于本文所提出的算法。文獻[17]提出雙線性池化模塊對網絡輸出特征進行融合，但由于自身主干網絡特征提取能力有限，而與本文算法性能差異最大。文獻[24]提出一種利用局部信息和全局外觀信息的非對稱多分支結構并引入弱監督增強學習能力，但是無法將類間特征分離且無弱監督能力，雖然有較大進步但實驗結果相比于本文算法下降了0.028 0。文獻[25]提出漸進訓練方法引導網絡學習，但特征融合時層間特征差距過大且部分特征未經充分提取也使得最終精度受到影響。文獻[26]采用多主干模型增強分類效果，最終結果依然弱于本文所提算法。文獻[27]融合多種FGC策略并通過度量學習獲得更具鑒別力的特征以幫助最終分類。對比上述研究，由分類結果可知，本文方法在IDRiD數據集中取得更有競爭力的性能結果。5折交叉驗證（交叉驗證-1～交叉驗證-5）訓練圖及學習率變化趨勢，如圖6所示。
圖6
學習率變化圖及模型交叉驗證QWK值變化趨勢圖
Figure6.
Learning rate change and model cross-validation QWK value change trend
如表3所示，現有FGC算法和最新的DR分級算法在Messidor-2數據集中，本文Sen對比文獻[28]和文獻[29]分別提高0.108 0和0.202 3。但是特異性低于文獻[29]，AUC均高于其他兩種算法，分別增長0.032 4和0.006 7。這也表明本文算法在Messidor-2數據集上依然具有一定的泛化能力。
4 結束語
本文提出了一個新的CHBP的DR-FGC模型，其整體結構由主干網絡（SEResNeXt）和CHBP分類器組成，并引入CL函數和FL函數緩解數據不均衡和種類差異度小的問題，同時利用漸進訓練方法提升精度。實驗結果證明，CHBP算法能夠有效地根據DR不同種類之間微小差異判斷DR等級。由于本文網絡整體參數量大，部分網絡還有精減的空間，而且引入漸進訓練方法大幅度增加了訓練時間，因此下一步研究重心將集中在精簡網絡和加快推理上。本文模型不但提高了DR分級的效率，而且不會出現醫生由于主觀原因而導致誤判，在DR分級領域具有很大的前景。
重要聲明
利益沖突聲明：本文全體作者均聲明不存在利益沖突。
作者貢獻聲明：梁禮明主要負責項目主持、平臺搭建、算法程序設計、協調溝通以及計劃安排，彭仁杰主要負責實驗流程、數據記錄與分析、論文編寫以及算法程序設計；馮駿和尹江主要負責提供實驗指導，數據分析指導，論文審閱修訂。
