陳卓 鄧興廣 繆曉強 鄭晗

1 廣東省藥品監(jiān)督管理局審評認證中心 （廣東 廣州 510080）

2 深圳技術大學4D醫(yī)學動態(tài)與功能成像聯(lián)合實驗室 （廣東 深圳 518055）

3 深圳市安健科技股份有限公司 （廣東 深圳 518055）

內容提要：多功能數(shù)字化動態(tài)數(shù)字X射線攝影（Digital Radiography，DR）具有高清影像、適時點片、實時回放、連續(xù)點片等多種功能，在臨床診斷中得到廣泛應用。動態(tài)DR可獲得肺部呼吸運動變化的圖像，對呼吸康復的效果測定、肺癌、慢性阻塞性肺疾病、支氣管哮喘、間質性肺炎等肺部疾病的診斷起到有用的價值，還可以觀察到肺部的血流分布和通氣分布情況。文章基于動態(tài)DR全呼吸周期的圖像，通過傳統(tǒng)方法和深度學習的結合進行量化分析，包括肺野識別、左右肺劃分、肺部面積計算、膈肌定位、實時輸出肺面積變化曲線和膈肌運動曲線。為臨床醫(yī)生提供多項量化指標，并且為臨床動態(tài)DR診斷提供輔助支持。

隨著生活質量日益提高，人們更加重視自身的健康狀況，對醫(yī)療資源的需求也逐漸增大。目前，大多數(shù)人體檢采用的靜態(tài)數(shù)字X射線攝影（Digital Radiography，DR）拍攝，動態(tài)DR拍攝尚未普及。動態(tài)DR有著巨大的潛力，通過對動態(tài)DR影像的分析可以觀察到肺部血流和通氣的分布情況，以及肺面積、橫膈膜和肺密度的動態(tài)變化情況。本文分別采用傳統(tǒng)方法和深度學習方法研究動態(tài)DR圖像的肺分割以及左右橫膈膜的獲取，目標是達成對動態(tài)DR影像肺區(qū)域的實時分割以及得到左右肺面積和橫膈膜位置的實時變化曲線。

1.國內外研究進展

1.1 肺分割

基于動態(tài)DR的肺分割，可通過單幀的肺分割加上幀與幀之間的平滑實現(xiàn)。為了精確地分割肺野，在過去的四十年里，人們已經(jīng)發(fā)展出了許多方法[1]。根據(jù)Mittal等[1]的研究，他們對肺分割的方法進行了全面的調查，重點關注了它們的基本原理、所使用的數(shù)據(jù)集、分割的性能和相對的優(yōu)缺點。Van Ginneken等[2]將肺野分割方法分為四大類：①基于規(guī)則的方法；②基于像素分類（PC）的方法；③基于可變形模型的方法；④混合方法。Van Ginneken等[2]調查的目的是對胸部圖像計算機分割的文獻進行分類和簡要回顧，其中包括了過去30年發(fā)表的150多篇論文。

1.2 橫膈膜運動分析

Baltruschat等[3]通過對172名志愿者的實驗，使用動態(tài)DR對站立呼吸時膈肌的時間分辨率定量分析，結果表明，在站立式潮汐呼吸時，橫膈肌的平均移動距離分別為11.0mm（右）和14.9mm（左）。左側膈肌運動明顯大于右側。較高的BMI和潮氣量與雙側橫膈肌的短途活動增加有關。本文將基于此項研究繼續(xù)探索橫膈膜的運動問題。

2.動態(tài)DR圖像的傳統(tǒng)方法肺分割

本文圍繞動態(tài)DR圖像的肺分割，分別采用傳統(tǒng)方法進行處理。包括基于閾值的肺分割的全閾值分割、最大類間方差法分割、自適應閾值法分割和基于非剛性配準的分割。

2.1 基于傳統(tǒng)方法的肺分割

圖像閾值分割是一種傳統(tǒng)的最常用的圖像分割方法，因其實現(xiàn)簡單、計算量小、性能較穩(wěn)定而成為圖像分割中最基本和應用最廣泛的分割技術。它特別適用于目標和背景占據(jù)不同灰度級范圍的圖像。它不僅可以大大壓縮數(shù)據(jù)量，而且也大大簡化了分析和處理步驟，因此，在很多情況下圖像閾值分割是進行圖像分析、特征提取與模式識別之前的必要的圖像預處理過程。

圖像閾值的目的是要按照灰度級，對像素集合進行一個劃分，得到的每個子集形成一個與現(xiàn)實景物相對應的區(qū)域，各個區(qū)域內部具有一致的屬性，而相鄰區(qū)域不具有這種一致屬性。這樣的劃分可以通過從灰度級出發(fā)選取一個或多個閾值來實現(xiàn)，從而實現(xiàn)目標與背景的分離。在醫(yī)學圖像上而言，組織器官的病灶的分割給臨床診斷提供了一種有效的輔助手段，下文將分別介紹三種常見的基于閾值的肺分割算法。

2.1.1 全閾值分割

全閾值分割指將灰度值大于閾值的像素設為一種顏色，小于或等于閾值的像素設為另外一種顏色，在OpenCV庫中可以實現(xiàn)全閾值分割。

用數(shù)學表達式來表示，則可設T為閾值，使得原始圖像f(x,y)轉換為二值圖g(x,y)，其中x和y表示像素點的位置，全閾值分割圖像時則滿足公式（1）。

基于全閾值肺分割的總體方案如圖1所示，并且每個步驟的結果如圖2所示，并且制作了Matlab可操作交互式UI界面，界面中包含左右肺面積變化曲線的實時分割結果，UI界面如圖3所示。

圖1. 基于全閾值肺分割的總體方案

圖2. 基于閾值的肺分割的分步結果

圖3. 肺分割的Matlab UI界面

2.1.2 最大類間方差法分割

最大類間方差法（OTSU）也叫大津法，是一種確定圖像二值化分割閾值的算法，由日本學者大津于1979年提出。從大津法的原理上來講，該方法又稱作最大類間方差法，因為按照大津法求得的閾值進行圖像二值化分割后，前景與背景圖像的類間方差最大。

OTSU是一種基于全局的二值化算法，它是根據(jù)圖像的灰度特性，將圖像分為前景和背景兩個部分。當取最佳閾值時，兩部分之間的差別應該是最大的，在OTSU算法中所采用的衡量差別的標準就是較為常見的最大類間方差。前景和背景之間的類間方差如果越大，就說明構成圖像的兩個部分之間的差別越大，當部分目標被錯分為背景或部分背景被錯分為目標，都會導致兩部分差別變小，當所取閾值的分割使類間方差最大時就意味著錯分概率最小。

記T為前景與背景的分割閾值，前景點數(shù)占圖像比例為w0，平均灰度為u0；背景點數(shù)占圖像比例為w1，平均灰度為u1，圖像的總平均灰度為u，前景和背景圖像方差為g，具體計算見公式（2）～（4）。

由公式（2）和（3）可得公式（4）：

當方差g最大時，可以認為此時前景和背景差異最大，此時的灰度T是最佳閾值。類間方差法對噪聲以及目標大小十分敏感，它僅對類間方差函數(shù)為單峰的圖像產(chǎn)生較好的分割效果。當目標與背景的大小比例懸殊時（例如受光照不均、反光或背景復雜等因素影響），類間方差函數(shù)可能呈現(xiàn)雙峰或多峰，此時效果不好。

在OpenCV中實現(xiàn)OTSU算法使用的是cv.threshold()函數(shù)，指定使用cv.THRESH_OTSU即可。

2.1.3 自適應閾值法分割

自適應閾值分割也叫做局部閾值化，它是根據(jù)像素的鄰域塊的像素值分布來確定該像素位置上的閾值。這樣做的好處在于每個像素位置處的閾值不是固定不變的，而是由其周圍鄰域像素的分布來決定的。亮度較高的圖像區(qū)域的閾值通常會較高，而亮度較低的圖像區(qū)域的閾值則會相適應地變小。不同亮度、對比度、紋理的局部圖像區(qū)域將會擁有相對應的局部閾值。

常用的局部自適應閾值有局部鄰域塊的均值和局部鄰域塊的高斯加權和。一種較為簡單的自適應閾值選取方法是對每個像素確定以其自身為中心的一個窗口，尋找窗口內像素的最大值與最小值，并取二者的平均值作為閾值，或者將窗口內所有像素的平均值作為閾值，亦或者將窗口內的所有像素的高斯卷積作為閾值。

全局閾值化使用唯一的閾值，對于亮度分布差異較大的圖像，難以找到一個合適的閾值。而使用自適應的閾值分割時，因為閾值是自適應的，所以可以將物體分割出來，效果要優(yōu)于全閾值分割。

2.2 基于非剛性配準的肺分割

此方法是Sema Candemir等[4]提出來的一種巧妙的方法，使用到了非剛性配準技術，和預先帶有人工分割標簽的數(shù)據(jù)集。非剛性配準系統(tǒng)由以下三個階段組成。

2.2.1 Ⅰ階段

在訓練數(shù)據(jù)庫中識別5個與患者最相似的圖像，并使用這個訓練圖像和訓練圖像的標簽來扭曲出患者特定的肺邊界。通過比較患者X射線與人工分割的肺圖像，能顯著改善之后圖像配準的準確率。

首先使用Partial Radon Transforms，它是一個積分變換，將定義在二維平面上的一個函數(shù)f(x,y)沿著平面上的任意一條直線做線積分，來比較和排序兩個患者的肺圖像之間的相似性，Partial Radon Transforms如公式（5）所示。

注：I(·)是二維圖像，δ(·)是二維脈沖函數(shù)，ρ=xcosθ+ysinθ，θ是法向量和x軸的夾角。

訓練數(shù)據(jù)庫中的所有圖像預先計算水平和垂直投影輪廓。計算直方圖在垂直和水平方向上的投影，然后用Bhattacharyya系數(shù)來測量訓練數(shù)據(jù)庫和患者胸部X射線之間的每個投影輪廓的相似性，計算相似性如公式（6）所示。

注：p1(x)、p2(x)是圖像I1和I2的水平投影，q1(x)、q2(x)是圖像I1和I2的垂直投影，x和y分別是投影配置文件的直方條，n和m是投影配置文件的直方條的數(shù)量，α=n/(n+m)是每個配置文件的相對比重。

2.2.2 Ⅱ階段

選用SIFT特征提取算法進行特征點提取。該算法實現(xiàn)的主要步驟包括：檢測尺度、計算空間極值點、確定特征點位置、提取特征點主方向、生成SIFT特征描述符。

在I階段從數(shù)據(jù)庫找到與患者最相似的數(shù)據(jù)，分別在數(shù)據(jù)庫和患者的X射線圖上找出SIFT特征點，通過對應匹配之間的空間距離移動所有像素，將變換映射應用于所有像素，在人工分割標簽的基礎上扭曲出患者特定的掩碼圖。圖像扭曲獲取掩碼圖的公式如（7）所示。

注：P是X射線的像素集，N是空間鄰域集（Spatial Neighborhood Set），S1、S2是由SIFT描述符向量表示的SIFT圖像，是關于P的光流向量（Flow Vectors），t和d是截斷的閾值。

2.2.3 Ⅲ階段

該系統(tǒng)利用圖像特性和上一階段計算出的肺模型概率，來檢測X射線圖像的肺邊界。使用圖割算法進行圖像分割，并用目標函數(shù)建模分割過程。

綜上所述，本文通過借鑒Sema Candemir等[4]使用非剛性配準的數(shù)據(jù)集在胸片中進行肺分割，將此技術運用到動態(tài)DR圖像中的肺分割中，完成了現(xiàn)有十四套DR數(shù)據(jù)肺區(qū)域的實時分割。但是由于數(shù)據(jù)集較小，魯棒性不是特別好，并且分割速度較慢。

部分結果的UI界面如圖4所示，界面中有四個坐標圖，左上表為左肺面積變化曲線、右上表為右肺面積變化曲線、左下角為左橫膈膜位移變化曲線、右下角為右橫膈膜位移變化曲線。

圖4. 基于非剛性配準的肺分割結果

根據(jù)分割的結果來看，部分數(shù)據(jù)出現(xiàn)突變點或者分割效果不穩(wěn)定，這與標簽的數(shù)量有關，標簽的多樣性和豐富性對結果有影響，后期可以通過添加數(shù)據(jù)集中的樣本數(shù)量和多樣化來改善分割效果。

基于非剛性配準肺分割的結果，可以提取到橫膈膜的位移信息。以第一幀圖像橫膈膜的最高點為基準，通過固定第一幀的x軸，來確定后續(xù)幀的橫膈膜位置，位置結果如圖5所示。

圖5. 橫膈膜位置的自動選取

3.基于深度學習的肺分割

3.1 環(huán)境搭建

本文采用基于U-Net改編的網(wǎng)絡進行實驗，具體實驗平臺的相關配置如表1所示。

表1. 實驗平臺的相關配置

PyTorch是一個開源的Python機器學習庫，用于自然語言處理、深度學習等領域。2017年1月，由Facebook人工智能研究院基于Torch推出了PyTorch。它是一個基于Python的可續(xù)計算包，具有強大的GPU加速的張量計算（如NumPy），并且包含自動求導系統(tǒng)的深度神經(jīng)網(wǎng)絡。

計算統(tǒng)一設備體系結構（Compute Unified Device Architecture，CUDA），是顯卡廠商NVIDIA推出的運算平臺。CUDA是一種通用并行計算架構，該架構使GPU能夠解決復雜的計算問題。

cuDNN是用于深度神經(jīng)網(wǎng)絡的GPU加速庫。它強調性能、易用性和低內存開銷。NVIDIA cuDNN可以集成到更高級別的機器學習框架中，如Tensor flow和Caffe。簡單的插入式設計可以讓開發(fā)人員專注于設計和實現(xiàn)神經(jīng)網(wǎng)絡模型，而不是簡單調整性能，同時還可以在GPU上實現(xiàn)高性能現(xiàn)代并行計算。

3.2 數(shù)據(jù)準備

本文采用的數(shù)據(jù)集來自Kaggle網(wǎng)站上關于深圳醫(yī)院的肺部DR圖像，數(shù)據(jù)集包含已經(jīng)分割好的標簽，為了提高模型的魯棒性和泛化能力，本研究在現(xiàn)有的樣本上進行數(shù)據(jù)增強。

數(shù)據(jù)增強也叫數(shù)據(jù)擴增，即在不實質性的增加數(shù)據(jù)的情況下，讓有限的數(shù)據(jù)產(chǎn)生等價于更多數(shù)據(jù)的價值。數(shù)據(jù)增強可以分為有監(jiān)督的數(shù)據(jù)增強和無監(jiān)督的數(shù)據(jù)增強方法。其中，有監(jiān)督的數(shù)據(jù)增強又可以分為單樣本數(shù)據(jù)增強和多樣本數(shù)據(jù)增強方法，無監(jiān)督的數(shù)據(jù)增強分為生成新的數(shù)據(jù)和學習增強策略兩個方向。

本文采用的數(shù)據(jù)增強的方法包括隨機旋轉、隨機裁剪和均衡直方圖化。通過調用PyTorch框架里的torchvision.transforms來實現(xiàn)，然后對數(shù)據(jù)增強后的圖像進行訓練。

數(shù)據(jù)增強的本質是為了增強模型的泛化能力，它與其他的一些方法比如DropOut，權重衰減有一定的區(qū)別。①權重衰減、DropOut等方法是專門設計來限制模型的有效容量的，用于減少過擬合，這一類是顯式的正則化方法。研究表明，這一類方法可以提高泛化能力，但并非必要且能力有限，同時參數(shù)高度依賴于網(wǎng)絡結構。②數(shù)據(jù)增強，沒有降低網(wǎng)絡的容量，也不增加計算復雜度和調參工程量，是隱式的規(guī)整化方法，在實際應用中更具有意義。

3.3 BD-U-Net網(wǎng)絡結構

本文模型基于U-Net框架在卷積塊中加入了批標準化（Batch Normalization）和DropOut來提高模型的泛化能力。U-Net模型是菲茲堡大學Olaf Ronneberger等[5]在2015年參加ISBI競賽時，提出的一種全卷積神經(jīng)網(wǎng)絡，最初用于完成細胞分割任務。U-Net模型結構優(yōu)雅精巧，對于數(shù)據(jù)較少的醫(yī)學圖像，端到端的分割任務效果較好。

Olaf Ronneberger等[5]提出的U-Net的U形結構如圖2～6所示。網(wǎng)絡是一個經(jīng)典的全卷積網(wǎng)絡（即網(wǎng)絡中沒有全連接操作）。網(wǎng)絡的輸入是一張572×572的邊緣經(jīng)過鏡像操作的圖片。

圖6. U-Net網(wǎng)絡結構[10]

網(wǎng)絡的左側在論文中將這一部分叫做壓縮路徑。壓縮路徑由4個塊組成，每個塊使用了2個有效卷積和1個Max Pooling降采樣，每次降采樣之后Feature Map的個數(shù)乘2，因此有了圖中所示的Feature Map尺寸變化。最終得到了尺寸為32×32的Feature Map。

網(wǎng)絡的右側部分在論文中叫做擴展路徑。同樣由4個塊組成，每個塊開始之前通過反卷積將Feature Map的尺寸乘2，同時將其個數(shù)減半（最后一層略有不同），然后和左側對稱的壓縮路徑的Feature Map合并，由于左側壓縮路徑和右側擴展路徑的Feature Map的尺寸不一樣，U-Net是通過將壓縮路徑的Feature Map裁剪到和擴展路徑相同尺寸的Feature Map進行歸一化的（即圖2～6中左側虛線部分）。擴展路徑的卷積操作依舊使用的是有效卷積操作，最終得到的Feature Map的尺寸是388×388。由于該任務是一個二分類任務，所以網(wǎng)絡有兩個輸出Feature Map。

本文在卷積過后增加Batch Normalization和dropout操作，以此提高模型的泛化能力，并且防止模型過擬合。其中為了簡潔表格，定義了兩個模塊，分別為Conv_BN_DP和ConvT_BN_DP，如表2所示。修改過的模型網(wǎng)絡結構及每層的參數(shù)和每層輸出圖像尺寸如表3所示，并且BD-U-Net網(wǎng)絡的總計算量如表4所示。

表2. BD-U-Net模塊

表3. BD-U-Net 結構

表4. BD-U-Net總參數(shù)

Batch Normalization是Sergey Ioffe等[6]在2015年研發(fā)的，普通的歸一化操作（Z-score標準化）歸一化到均值為0、方差為1的公式如下，即每個樣本減去樣本數(shù)據(jù)的平均值后再除以樣本數(shù)據(jù)的標準差，如公式（8）所示。

注：ε是個很小的數(shù)，防止分母為0，γ和β是可訓練的參數(shù)。

3.4 訓練與驗證

本文采用交叉熵損失函數(shù)，Golik等[7]提出在隨機初始化權重的情況下，基于平方誤差的神經(jīng)網(wǎng)絡不能收斂到一個好的局部最優(yōu)解，而交叉熵更加快速地使網(wǎng)絡收斂。

學習率是一個至關重要的參數(shù)，影響模型的收斂速度，本文使用動態(tài)學習率，通過調用torch.optim.lrscheduler.StepLR，可以實現(xiàn)每30輪訓練學習率遞減為原來的10%，本文采用的訓練參數(shù)如表5所示。

表5. 訓練參數(shù)

本研究制作的數(shù)據(jù)集包含700張肺部圖片，將其分為訓練集：驗證集：測試集=7:2:1。每訓練完一輪后，就進行一次模型效果的驗證，全部輪數(shù)訓練完后，再用測試集中的數(shù)據(jù)進行測試。訓練損失與訓練輪數(shù)的關系如圖7所示，驗證準確率與訓練輪數(shù)的關系如圖8所示。

圖7. 訓練損失與訓練輪數(shù)的關系圖

圖8. 驗證準確率與訓練輪數(shù)的關系

3.5 測試結果

將測試數(shù)據(jù)輸入模型后得到肺的掩碼圖，通過圖像形態(tài)學腐蝕、膨脹的操作，可以光滑肺邊緣和去除小連通域，也可以通過掩碼圖減去膨脹的掩碼圖，得到肺掩碼圖的輪廓線。肺部分割效果如圖9所示。

圖9. 肺部分割效果

通過計算每一幀肺部圖像的肺面積和肺密度，就可以的到肺面積變化曲線和肺密度變化曲線，曲線的形態(tài)以及動態(tài)區(qū)間，可以輔助醫(yī)生進行診斷。動態(tài)曲線如圖10所示。

圖10. 肺部動態(tài)曲線

綜上所述，通過傳統(tǒng)方法和深度學習的肺分割結果來看，深度學習的速度比傳統(tǒng)方法快，其魯棒性比傳統(tǒng)方法好，其準確率比傳統(tǒng)方法高。

4.小結

醫(yī)學影像是輔助臨床疾病診療重要手段，其中DR是最為普及的醫(yī)療儀器之一，動態(tài)DR在此基礎上額外提供了不同時間點的動態(tài)信息。深度學習的模型對復雜圖像有高效的處理能力，在肺區(qū)域分割過程中，能夠在多幀DR處理時提高效率，并且準確率和魯棒性都很好。本研究實現(xiàn)了基于動態(tài)DR影像的實時左右肺分割以及橫膈膜提取，輸出了動態(tài)肺面積曲線，以及橫膈膜運動曲線。此曲線能夠輔助診斷患者的肺容量，可以輔助排查出慢性阻塞性肺疾病等肺部疾病，提升臨床診療質量以及效率。