基于LUT的快速3D氣道樹骨架線提取

劉明威， 顧力栩

(上海交通大學(xué) 生物醫(yī)學(xué)工程學(xué)院，上海 200240)

基于LUT的快速3D氣道樹骨架線提取

劉明威， 顧力栩

(上海交通大學(xué) 生物醫(yī)學(xué)工程學(xué)院，上海 200240)

為提高肺部支氣管骨架線的提取效率，提出并使用一種基于look-up-table(LUT) 的腐蝕細化算法. 分析建立腐蝕模型，并根據(jù)該模型優(yōu)化建立了LUT，以該LUT為依據(jù)通過索引查找對原始數(shù)據(jù)進行快速腐蝕細化，對得到的腐蝕結(jié)果進行剪枝處理以得到最終的骨架線. 實驗結(jié)果表明：提取過程中LUT的應(yīng)用從根本上降低了腐蝕細化中判斷的復(fù)雜度，將復(fù)雜的簡單點判斷問題轉(zhuǎn)化為LUT中的查詢問題，從而極大地優(yōu)化了腐蝕細化中關(guān)鍵的腐蝕過程. 相比傳統(tǒng)方法，基于LUT的腐蝕細化算法顯著提高了骨架線的提取速度，較傳統(tǒng)細化法提速近22.95倍.

計算機輔助診斷； 氣管樹； 虛擬支氣管鏡；骨架線；查找表

肺癌是威脅人類生命健康的惡性腫瘤疾病，中國目前總體肺癌發(fā)病率增長了19.48%，這很可能是加劇的環(huán)境污染問題如大規(guī)模霧霾引起的[1]. 早期診斷能大幅度提高病患存活率，但目前早期診斷率僅15%[2]. 常見肺癌早期診斷手段包括CT成像與光學(xué)內(nèi)窺鏡. CT成像擁有密度分辨率高、操作簡單等優(yōu)勢，但傳統(tǒng)CT只能提供單層的斷面信息，不直觀，且在診斷過程中病人需承受大量放射輻射；光學(xué)內(nèi)窺鏡可在零輻射下觀察到管腔內(nèi)部真實紋理色彩，但操作上要比CT成像復(fù)雜得多.

隨著計算機輔助診斷技術(shù)的發(fā)展，虛擬支氣管鏡VB[3](virtual bronchoscopy)技術(shù)在肺癌診斷領(lǐng)域正逐漸得以普遍應(yīng)用. 該技術(shù)利用從CT數(shù)據(jù)中重建的氣管模型，對肺部疾病診斷有重大意義. 臨床表明VB不僅可以顯示有無氣道異常如支氣管狹窄等，還可以配合導(dǎo)航系統(tǒng)引導(dǎo)經(jīng)支氣管鏡肺活檢，明顯提高支氣管鏡在肺部外周病灶診斷效率. 相對于傳統(tǒng)內(nèi)窺鏡而言，VB引導(dǎo)下的導(dǎo)航技術(shù)具有以下優(yōu)點：1)定位準確，診斷效能高；2)可快速引導(dǎo)支氣管鏡，縮短檢查時間，減少病患痛苦；3)可減少支氣管鏡檢查次數(shù)；4)可提供自由的虛擬觀測角度，降低支氣管異常的辨別難度； 5)可部分代替X線引導(dǎo)的支氣管鏡活檢，減少輻射；6)可避免盲目穿刺，降低并發(fā)癥發(fā)生率. 大量文獻證實虛擬內(nèi)窺鏡能在不失診斷準確率的基礎(chǔ)上，保障診斷的微創(chuàng)性[4-5].

在虛擬內(nèi)窺鏡實現(xiàn)過程中，快速的骨架線提取一直是其中的難點和要點. 而在VB中，一個良好提取的骨架線不僅可以提供直觀的支氣管拓撲結(jié)構(gòu)，還可以為管腔內(nèi)瀏覽起到引導(dǎo)作用. 與人體其他管腔如結(jié)腸、血管相比，支氣管結(jié)構(gòu)更為復(fù)雜[6]，體現(xiàn)在：1)復(fù)雜的拓撲結(jié)構(gòu)，擁有多級分支，多通路；2)各級分支粗細不一；3)個體差異性，不同人支氣管雖在級數(shù)上有一定相似性，但大小和分支位置卻存在差異. 目前三維骨架提取方式主要有距離變換法、中軸變換法和腐蝕細化算法. 然而由于支氣管的復(fù)雜性，這些方式在支氣管骨架線的提取上都略顯性能不足. 隨著當(dāng)今4DCT[7]數(shù)據(jù)的逐步普及，傳統(tǒng)方法的運算速度更顯得不盡人意.

本文提出一種基于腐蝕細化的改進三維骨架線提取方法，以克服傳統(tǒng)腐蝕細化方法中復(fù)雜度高、耗時長的問題. 本方法充分利用了當(dāng)今計算機多核以及高內(nèi)存的優(yōu)勢，在傳統(tǒng)三維腐蝕過程中引入LUT. 相關(guān)實驗證實，本方法能夠顯著提高骨架線的提取速度.

1 骨架線基本定義與現(xiàn)有方法綜述

1.1 骨架線基本定義

骨架(Skeleton)是指與原物體具有一致連通性和拓撲結(jié)構(gòu)的細曲線的一種理想表示，就是位于物體內(nèi)部，且能體現(xiàn)其拓撲特征的簡化圖形[8]. 在三維情況下提取的骨架線需要滿足如下要求[9]：1) 連通性，骨架線應(yīng)與目標(biāo)物體擁有相同的連通性；2) 居中性 ，骨架線盡量保持在管腔的中央；3) 單體素性，組成骨架線的體素要求為單元寬度. 對于簡單的柱狀組織，中心線只有唯一的一條，而對于支氣管樹或血管樹這種復(fù)雜的管狀結(jié)構(gòu)，中心線會有多條，也被稱作為骨架線.

1.2 現(xiàn)有方法綜述

就目前而言，二維的中心線提取算法已經(jīng)相對成熟. 相對而言，三維中心線的提取方式只有距離變換法、中軸變換法和腐蝕細化算法.

距離變換法[10]是指根據(jù)不同目的，將物體的體素點進行距離標(biāo)識. 目前主流的方法是基于雙距離場的方法[11]：首先對物體進行邊界距離變換，得到體素到達邊界的邊界距離場DFB (distance from boundary)；其次指定源點，進行源點距離變換，得到源距離場DFS (distance from source)；最后通過邊界距離場中局部最大值的篩選，配合源距離場進行骨架點的連接. 中軸變換法[12]在二維上可以用最大內(nèi)切圓盤的相關(guān)術(shù)語定義，而在三維上則表現(xiàn)為最大內(nèi)切球中心的集合. 其中圓盤中心到達邊界的距離稱為圓盤法線，一個目標(biāo)內(nèi)點的最大內(nèi)切圓在邊界上至少有兩個切點，而每個對稱點有兩個或多個圓盤法線. 腐蝕細化算法[13]的思想則更為樸素，根據(jù)物體歐拉特性[14]以及連通性的不變性，給出在刪除后不影響原物體拓撲結(jié)構(gòu)的簡單點，并以此為依據(jù)均勻、對稱地從不同方向?qū)δ繕?biāo)物體邊緣上的簡單點(Simple Voxel, SV)[15]進行剝離，直至僅剩下單像素寬的體素集合，去除該集合中的偽骨架點后得到最終的骨架線.

以上方法各自均有不同的局限性和一些共同不足. 距離變換法中，在構(gòu)建雙距離場時不僅因需要對各體素點以迭代的方式進行距離標(biāo)記而導(dǎo)致高時耗，還因需要與目標(biāo)物體同樣大小的內(nèi)存空間儲存距離場而帶來高內(nèi)存消耗. 同時，依據(jù)不同的歐氏距離定義建立的距離場會產(chǎn)生不同的候選中心點，直接影響骨架線結(jié)果的一致性[16]，且候選中心點通常是不連續(xù)的，連接時會造成骨架線在曲率較大位置的離心偏移，難以保證居中性. 中軸變換法中，由于中軸變換是理想的數(shù)學(xué)模型描述，主要基于數(shù)學(xué)的幾何推理證明，從而導(dǎo)致使用機器語言進行實現(xiàn)的難度大、算法復(fù)雜度高，常見的模型有地表火模型[17]，距離曲面脊線模型[18]. 以上兩種方法存在共同的缺陷就是易受邊緣噪聲干擾，且難以保證骨架的連續(xù)性. 腐蝕細化算法的不足在于需要不斷地對邊緣的簡單點進行判斷，現(xiàn)有的判斷模型復(fù)雜度均不夠理想，導(dǎo)致當(dāng)目標(biāo)物體體素點過多時，帶來巨大的時耗.

2 改進的腐蝕細化方法

考慮到基于距離場的算法對于邊緣噪聲的敏感性和內(nèi)存的大量消耗，以及傳統(tǒng)腐蝕細化算法中對邊緣點反復(fù)判斷的高復(fù)雜度，本文提出一種基于LUT的快速三維腐蝕細化算法. 該方法保留了傳統(tǒng)腐蝕細化算法的穩(wěn)定性，以及結(jié)果中對目標(biāo)物體拓撲結(jié)構(gòu)和連通性的完好保留性，并通過使用LUT的思想簡化細化過程，避免重復(fù)判斷，顯著提高了細化速度. 該方法首先從原始數(shù)據(jù)中用區(qū)域增長的方式分割出支氣管二值圖像，建立腐蝕細化模型并建立LUT，然后查找LUT進行腐蝕細化，最后對細化的結(jié)果進行剪枝處理.

2.1 支氣管二值數(shù)據(jù)的獲取

采用有監(jiān)控的動態(tài)閾值區(qū)域生長配合形態(tài)學(xué)閉操作來獲取組成支氣管，流程及效果如圖1所示. 首先以氣管壁與空氣CT值的差異為依據(jù)，通過應(yīng)用DMQ (dynamic marking QUEUE)的數(shù)據(jù)結(jié)構(gòu)進行增長，提高了提取的速度. 在頂層CT中選擇氣管入口作為增長原點壓入DMQ，目標(biāo)點離開DMQ時對其6鄰域進行增長判斷，若符合則壓入DMQ. 在判斷目標(biāo)點是否歸入增長區(qū)域時，采用了18鄰域平均值比較的方式，最大可能增加分割目標(biāo)點數(shù). 對每一輪新增的點與上一輪增長結(jié)果進行比較，以進行溢出監(jiān)控. 該方法得到一個26連通的初步結(jié)果，為平滑邊緣、填充管腔的空隙，在其基礎(chǔ)上使用公式IMG ·M = (IMG ? M) ⊕M進行形態(tài)學(xué)閉操作，用同一模版M分別對IMG進行膨脹腐蝕操作得到最終的二值分割結(jié)果.

2.2 細化方法

細化過程具體分為4步：1)腐蝕模型的建立； 2)基于腐蝕模型建立LUT(提高速度的關(guān)鍵步驟)；3)基于LUT進行快速細化； 4)對腐蝕結(jié)果進行剪枝處理.

圖1 支氣管分割流程

2.2.1 腐蝕模型的建立

在三維空間V3中，設(shè)中心點為c(x,y,z),對于任意點p(i,j,k)有如下的鄰域定義：

Nb6(p) = {p(i,j,k)∈V3| |i-x| + |j-y|+ |k-z| = 1},

Nb26(p) = {p(i,j,k)∈V3| |i-x|<2 ∩|j-y|<2∩|k-z| <2},

Nb18(p) = {p(i,j,k)∈V3| |i-x| + |j-y|+ |k-z| < 3}∩ Nb26(p).

將屬于目標(biāo)物體的體素定義為前景點Sfg，其余定義為背景點Sbg. 在前景點集合中存在子集Ssub-fg，若其中所有元素均符合以下四種判斷準則[19]，則稱之為簡單點SV(simplevoxel).

a) ?Sfgin Nb26and Sum(Sfgin Nb26)>1；

b) ?Sbgin Nb6；

c) ?Sfgin Nb26and is 26Connected(Sfg)；

d) ?Sbgin Nb6and is 6Connected By Nb18(Sbg).

a、b條件可以直接通過相應(yīng)鄰域中的前景點、背景點分布獲得，而c、d條件的判斷才是關(guān)鍵.

對于條件c，提出一種26鄰域區(qū)域增長法，如圖2所示. 根據(jù)條件c的描述易證明：若Nb26中的前景點可以彼此26連通，那么在該Nb26中，一定有且只有一個26連通的前景點區(qū)域. 因此，取26鄰域中的任意Sfg，在該Nb26中進行26鄰域的前景點區(qū)域增長，得到增長結(jié)果集合S′. 將其與26鄰域中的前景點構(gòu)成的集合S進行比較，若兩者點數(shù)相同，則符合條件c.

在條件d的判斷上，并不能完全與條件c等同，因為在Nb18中的Sbg中，6連通的數(shù)量并不唯一. 充分考慮到這一點，利用Nb18中Sbg的6鄰域局部區(qū)域增長法，使用列表套列表的結(jié)構(gòu)來儲存連通的6鄰域數(shù)組. 在遍歷結(jié)束后，將Nb6中的背景點集合S與連通的6鄰域列表中的Si依次進行比較，得到是否可以刪除.

圖2 條件c判斷

2.2.2 鄰域分布分析與LUT的建立

分析腐蝕模型可以看出，判斷一個前景點是否為簡單點是由a、b、c、d這4個條件的滿足情況決定的. 進一步分析，條件a為非孤立點的判斷，取決于Nb26；條件b為邊界點的判斷，取決于Nb6；條件c為Nb26的背景點的26連通情況，取決于Nb26；條件d為Nb26前景點的18連通情況，同樣取決于Nb26. 因此可得出結(jié)論，對于任意體素點p，可根據(jù)其Nb26(p)分布情況VC(voxelcombination)得到其簡單性.

在二值圖像前提下，任意體素點只有0或1兩種取值. 如圖3所示，若對26鄰域指定排列順序，則可等同于一個26bits的二進制數(shù)，進而任意一種VC將對應(yīng)唯一的26bits二進制數(shù)，而相對的任意一個26bits二進制數(shù)也能唯一映射為一種VC. 基于以上分析，可聯(lián)想到建立LUT輔助腐蝕細化. 在腐蝕細化過程中，將目標(biāo)點的26鄰域VC轉(zhuǎn)化為26bits的索引，在LUT中進行查詢判斷. 文獻[20]顯示LUT在二維的細化中得到了很好的應(yīng)用，值得一提的是，由于三維計算量巨大，空間分布復(fù)雜性以及受拘束于過去計算機的內(nèi)存，三維LUT鮮有人嘗試.

為避免建立LUT過程中的重復(fù)判斷，對體素的鄰域的空間分布狀態(tài)進行特征分析. 對于任意立方體，給出如下對于相似體素分布SVC(similarvoxelcombination)的定義：

1) 通過任意旋轉(zhuǎn)變換可以完全重合；

2) 通過左右手坐標(biāo)系對換可以完全重合，即

圖3 26鄰域轉(zhuǎn)換

(1)

(2)

(3)

條件1)可解釋為依次以立方體每一面朝向的坐標(biāo)系方向為觀察方向，而后以該方向為軸向，使用如式(1)～(3)中對應(yīng)方向的旋轉(zhuǎn)矩陣對立方體進行圍繞該軸向的三維旋轉(zhuǎn)，為保證立方體空間位置不變，這里旋轉(zhuǎn)的角度應(yīng)該分別取0°、90°、180°、270°，這樣通過旋轉(zhuǎn)可得4*6 共24種SVC；條件2)可解釋為條件1)基礎(chǔ)上的鏡像變換，因此同樣存在24種SVC. 因此可得出結(jié)論，對于各向異性的立方體，存在24+24共48種SVC，應(yīng)用在二值化的Nb26情況下，粗略估算在全空間2^26中存在2^26/48約為1398101種SVC.

(4)

而在二值化的Nb26情況下，并非所有分布都是各向異性，存在著大量的自反對稱，經(jīng)過式(4)的計算，其中i表示Nb26中前景點的個數(shù)，R(i) 表示Nb26中重復(fù)的組合個數(shù)，由于二值互補性，因此只需計算13次(0～12個前景點的情況)，外加13個前景點的情況. 針對i取值不同時需要具體列出對稱情況進行計算，但由于其復(fù)雜的旋轉(zhuǎn)對稱性，需要使用計算機輔助邏輯進行計算，最終得到不同前景點個數(shù)下SVC如表1所示，SVC一共有1426 144種，并在實際運算遍歷中加以檢驗.

表1 26鄰域中前景體素點與SVC關(guān)系表

基于以上分析，本文提出多線程并行配合雙布爾數(shù)組的算法建立LUT，數(shù)組大小為2^26，下標(biāo)代表該索引值所唯一對應(yīng)的VC. 其中一個數(shù)組為LUT，是后續(xù)腐蝕細化中所依賴的參照，其取值表示該VC下的中心點是否為簡單點；另一個數(shù)組為MT(marking table)，其取值用來表示該VC是否已被計算. 將所有2^26個索引等分為多個索引段，平均分配給多個子線程，配合多核CPU進行計算，以達到最高的計算效率.

所有子線程將共享LUT和MT，每個子線程將同時并行遍歷被分配的索引段. 在具體計算索引時，每一次簡單點判斷的成本很高，但應(yīng)用雙數(shù)組算法可有效避免重復(fù)計算. 對于一個索引值，首先在MT中查找判斷該索引值是否已被計算，若尚未計算，則將該索引值轉(zhuǎn)換為其對應(yīng)的VC，并根據(jù)前文的判斷條件和方法得出其是否為簡單點；而后通過正方體的旋轉(zhuǎn)對稱性，獲取該VC的所有其他47個SVC，并轉(zhuǎn)換為索引序列IA(index array)；最后在LUT中將IA中的索引值賦予與該VC相同的簡單性，同時在MT中將IA中的這些索引值標(biāo)記為已判斷. 雙數(shù)組的優(yōu)勢也在此，通過標(biāo)記MT避免了遍歷到SVC所對應(yīng)索引值時潛在的重復(fù)簡單性判斷，從而大大減少了判斷的次數(shù)，整體加速了LUT的建立過程. 為更為直觀，下面給出該流程框架的偽代碼.

函數(shù) 1: CreateLUT(創(chuàng)建LUT 的主函數(shù))

輸入: 需要的線程數(shù)目 threadNumber.

輸出: 計算好的lutArray.

for eachiinthreadNumberdo

subThreadsi=start_new_thread(MultiThreadCalculate(lutArray,mtArray,rangei) )

endfor

wait_sub_threads(subThreads)

write_to_bin_file(lutArray)

函數(shù)2:MultiThreadCalCulate(用于計算指定區(qū)段的子函數(shù))

輸入: 公共使用的lutArray，用于標(biāo)記是否計算過的mtArray，需要遍歷的區(qū)段range.

輸出: 計算指定區(qū)段range后的lutArray， 標(biāo)記指定區(qū)段range后的mtArray

foreachiinRange

ifnotmtArrayithen

isSimple=isSimplePoint(i)

similarIndexArray=getSimilarIndex(i)

foreachjinsimilarIndexArray

lutArrayj=isSimple

mtArrayj=true

endfor

endif

endfor

2.2.3 基于LUT進行迭代腐蝕細化

在建立好LUT的基礎(chǔ)上，傳統(tǒng)的腐蝕細化過程便得以優(yōu)化，其中對于SV的判斷可以直接依據(jù)LUT的結(jié)果，而不是重新對前文4個定義條件的判斷. 在進行細化前先將圖像中的前景點存入臨時數(shù)組PointListfg，以避免重復(fù)遍歷，提高效率. 而后進行如下操作：

1)在PointListfg中根據(jù)Nb26體素分布情況尋找邊緣點，得到PointListborder.

2)遍歷PointListborder中的點，將Nb26體素分布轉(zhuǎn)換為26bits索引值，根據(jù)LUT的判斷值進行簡單點的判斷刪除.

3)重復(fù)步驟1)、2)直至PointListborder中再也無法找到簡單點為止.

2.2.4 后續(xù)剪枝處理

實驗證實，通過上述腐蝕細化算法所得結(jié)果會在一些較粗的主氣管上出現(xiàn)一些額外的毛刺，這無疑會對后續(xù)的內(nèi)窺鏡路徑的規(guī)劃帶來不必要的分枝，需要進行剪枝處理.

根據(jù)毛刺分支長度短的特點，設(shè)計了一種尺寸自適應(yīng)的剪枝算法. 記當(dāng)前的前景點集合為S，首先根據(jù)骨架線的體素數(shù)設(shè)置適當(dāng)?shù)募糁﹂撝礣hreshold，設(shè)閾值百分比為th%，這里th取經(jīng)驗值1.7，因此經(jīng)計算得到 Threshold=S×th%. 接著遍歷細化后的骨架線，得到所有端點PointListend，以PointListend中的點為起點進行遍歷，直至分支點，若其長度小于Threshold則將該支刪除.

3 實驗結(jié)果

3.1 實驗環(huán)境與實驗數(shù)據(jù)

實驗數(shù)據(jù)采用從上海肺科醫(yī)院獲得的5組肺CT圖像序列，每組序列均包含有最少265層到最多468層CT斷層圖像，每層圖像的分辨率均是512*512 且灰度值在-2048到4196區(qū)間內(nèi). 用于實驗的PCCPU為Intel(R)Core(TM)i7-3610QM，RAM為6GB，OS為Win7 64-bit. 算法的實現(xiàn)基于VTKv6.1.0開源庫，配合QT在VS2008IDE下進行開發(fā)和實驗. 3.2 LUT的建立

為獲取最適合本實驗環(huán)境的線程數(shù)量，以得到最快的LUT建立速度，以線程數(shù)作為輸入?yún)?shù)進行實驗，得到如圖4所示的結(jié)果.

圖4 LUT建立線程時間函數(shù)

可見當(dāng)線程數(shù)少于5個時，增加線程可顯著減少用時，5個線程的處理過程表現(xiàn)最為優(yōu)秀，當(dāng)線程多于5個時，用時無明顯減少.

3.3 骨架線提取結(jié)果與比較

為比較基于LUT的方法與傳統(tǒng)腐蝕細化算法的效率，采用如下的對比試驗，使用相同的數(shù)據(jù)集進行重復(fù)試驗100次，表2列出了實驗結(jié)果.

在腐蝕細化階段，基于LUT的平均時間為482.74 ms，而基于傳統(tǒng)腐蝕細化算法的平均時間為11 077.7 ms，可見本文方法提高了22.95倍速度.

3.4 骨架線的提取效果

圖5為本文算法提取骨架線的結(jié)果與基于距離變換算法[21]的結(jié)果，其中圖5(a)、(b)、(c)三組結(jié)果中，上方為本文算法的結(jié)果，下方為距離變換算法的結(jié)果. 通過對比可以看出，基于距離變換的方法會在末端發(fā)生斷裂，不連續(xù)，在同一支氣管段上可能產(chǎn)生多分支，不方便路徑規(guī)劃，且操作中需要手工選定端點，既不方便也無法保證連續(xù)性；而本文方法在支氣管骨架線的提取上既可以保證骨架線在細末端的連續(xù)性，不發(fā)生斷裂，也無需手動，具備可重復(fù)性. 圖5中的三組數(shù)據(jù)顯示，該算法很好地保持了支氣管拓撲結(jié)果，同時也保障骨架線的中心性、連續(xù)性、一致性和可重復(fù)性，為內(nèi)窺鏡的漫游提供可能.

表2 實驗結(jié)果對比表

(a) 512*512*265 (b) 512*512*429 (c) 512*512*438

4 實驗討論

4.1 基于LUT算法的優(yōu)越性

本文同時還與基于距離場[21]算法得出的速度結(jié)果進行比較，由于使用不同的數(shù)據(jù)集，因此定義指標(biāo)ts(thining speed)來表示算法的效率作為比較依據(jù)，ts表示毫秒內(nèi)處理的立方點數(shù), 單位為103voxels·ms-1. 經(jīng)計算，文獻[21]的速度為220.16 ts，而本文的速度為45 253.39 ts，提高了205倍. 結(jié)合本文與以上兩種優(yōu)化過的傳統(tǒng)骨架線提取方式[21-22]的比較，可以充分看出基于LUT的算法能在保證骨架線拓撲結(jié)構(gòu)完整的基礎(chǔ)上，顯著提高運行速度.

總的來說，基于LUT的算法在本質(zhì)上優(yōu)化了傳統(tǒng)細化過程中的腐蝕細化環(huán)節(jié)，從而在整體上提高了骨架線提取的速度. 不僅如此，LUT的另一大優(yōu)勢是將LUT的建立與LUT的參照分離，這樣在簡單點的定義發(fā)生變化時，可以直接更新LUT的數(shù)據(jù)，而不需要更改后續(xù)整個程序.

4.2 參數(shù)與復(fù)雜度分析

基于LUT的腐蝕細化算法中，主要的輸入?yún)?shù)體現(xiàn)于兩點：一是LUT建立過程中線程數(shù)的選擇，二是應(yīng)用LUT進行腐蝕細化時，不同輸入物體的體素數(shù)量.

在LUT建立過程中，線程數(shù)是重要的輸入?yún)?shù)，并影響最終LUT建成的時耗. 實驗結(jié)果表明5個線程為最優(yōu)，分析其原因為：當(dāng)線程池大小為CPU數(shù)+1時，可達到最優(yōu)利用率，即當(dāng)執(zhí)行密集型任務(wù)時既使偶爾因頁缺失故障或者其他原因?qū)е聲和r，此額外線程也可確保CPU的時鐘周期不被浪費. 實驗PC 為四核CPU，因此5個線程可達到最優(yōu)化. 算法的空間復(fù)雜度為2*2^26*1 byte = 128 MB，用于緩存雙數(shù)組中的LUT和MT，在LUT建成結(jié)束后，得到的結(jié)果為一個64MB的二進制文件，可以直接作為查找腐蝕的依據(jù). 在時間上，對腐蝕細化的判斷條件c的判斷由于應(yīng)用了26鄰域內(nèi)區(qū)域增長法，只用一重循環(huán)，是一個線性過程，將該步驟的時間復(fù)雜度降低為O(n)；而對于條件d的判斷，在區(qū)域增長上的平均時間復(fù)雜度依舊為O(n)，而在列表比較上的平均時間復(fù)雜度為O(n2)，兩步之間為串行，因而整體的復(fù)雜度依舊為O(n2)，但較傳統(tǒng)方法相比依舊有所提升.

在腐蝕細化的過程中，由于LUT已經(jīng)確立，因此主要的輸入?yún)?shù)為物體的體素數(shù). 對3.3中的實驗結(jié)果進行線性擬合，得到t=0.007 2v-152.91，擬合優(yōu)度R2為0.988 5，其中t表示腐蝕細化的時耗，單位為ms，v表示物體的體素數(shù). 因此可得，腐蝕細化的時耗正比于輸入物體的體素數(shù). 算法的空間復(fù)雜度取決于物體的體素數(shù)，會隨著CT采集時分辨率的不同而有所變化，另外需考慮LUT占用的64MB；在時間復(fù)雜度上，查詢LUT接近于O(1)，對體素點進行擦除的復(fù)雜度也為O(1)，因此整體復(fù)雜度仍為O(1).

4.3 基于LUT算法的不足

與基于其他理論的骨架線提取方式如距離變換法相比，基于LUT的腐蝕細化法不足在于最終得到的骨架線是組成骨架線的離散點，缺少點與點之間拓撲關(guān)系的描述，為后續(xù)骨架線拓撲結(jié)構(gòu)分析帶來了一定的困難. 同時剪枝過程也有一定的局限性，無法精確區(qū)分細末較短的支氣管與偽分支. 而在LUT的使用中依舊可以作適當(dāng)優(yōu)化，如在查詢表的尺寸上，若是對內(nèi)存有苛刻的要求，可將LUT壓縮為8MB，以bit的形式來儲存每一個索引的索引值. 由于對bit的檢索需要額外的比較，其效率可能會有一定下降，這是一個棘手的問題.

受實驗數(shù)據(jù)源限制，未進行大規(guī)模驗證. 同時受實驗室硬件條件限制，缺少在其他配置相當(dāng)?shù)挠嬎銠C上進行重復(fù)驗證，也沒有獲取該算法在更差運行環(huán)境中的表現(xiàn). 更加充足的數(shù)據(jù)源和豐富的實驗平臺將能夠?qū)υ撍惴ㄓ幸粋€更全面的評估，這也是未來工作的目標(biāo).

5 結(jié) 語

基于LUT的腐蝕細化算法在傳統(tǒng)的三維腐蝕細化算法中引入LUT的思想，通過鄰域區(qū)域增長法優(yōu)化了簡單點腐蝕模型的判斷，并通過二值化26鄰域旋轉(zhuǎn)對稱性的分析進行快速的LUT建立，在腐蝕細化中使用LUT取代傳統(tǒng)算法中對簡單點的復(fù)雜重復(fù)判斷過程，從而顯著提高了骨架線的提取速度. 實驗結(jié)果表明：本文方法在速度上是傳統(tǒng)細化法方法[22]的22.95倍，傳統(tǒng)距離場方法[21]的近220倍，可以完好地保留肺部氣管樹的拓撲結(jié)構(gòu)，具備易操作性、連續(xù)性和可重復(fù)性，為內(nèi)窺鏡瀏覽導(dǎo)航提供可能.

[1] CHEN W, ZHENG R, ZENG H, et al. Annual report on status of cancer in China, 2011[J]. Chinese Journal of Cancer Research, 2015, 27(1): 2.

[2] 韋春暉. 肺癌早期診斷進展[J]. 臨床肺科雜志, 2010, 15(8):1136-1138.

WEI Chunhui. Progress in early diagnosis of lung cancer[J]. Journal of Clinical Pulmonary Medicine,2010,15(8):1136-1138.

[3] JUNG S Y, PAE S Y, CHUNG S M, et al. Three-dimensional CT with virtual bronchoscopy: a useful modality for bronchial foreign bodies in pediatric patients[J]. European Archives of Oto-Rhino-Laryngology, 2012, 269(1): 223-228.

[4] WERNER H, SANTOS J R L D, FONTES R, et al. Virtual bronchoscopy for evaluating cervical tumors of the fetus[J]. Ultrasound in Obstetrics & Gynecology, 2013, 41(1): 90-94.

[5] DE WEVER W, VANDECAVEYE V, LANCIOTTI S, et al. Multidetector CT-generated virtual bronchoscopy: an illustrated review of the potential clinical indications[J]. European Respiratory Journal, 2004, 23(5): 776-782.

[6] DENIZ A, HOFFMAN E A, GEOFFREY M L, et al. Segmentation and analysis of the human airway tree from three-dimensional X-ray CT images.[J]. IEEE Transactions on Medical Imaging, 2003, 22(8):940-50.

[7] COLE A J, O'Hare J M, McMahon S J, et al. Investigating the potential impact of four-dimensional computed tomography (4DCT) on toxicity, outcomes and dose escalation for radical lung cancer radiotherapy[J]. Clinical Oncology, 2014, 26(3): 142-150.

[8] 陳剛, 呂煊, 王志成, 等. 肺 CT 圖像的血管骨架化方法[J]. 計算機科學(xué), 2013, 40(5): 274-278.

CHEN Gang, LV Xuan, WANG Zhicheng, et al. Vessel Skeletonization Method for Lung CT Images[J].Computer Science, 2013, 40(5): 274-278.

[9] ZHOU Y, TOGA A W. Efficient skeletonization of volumetric objects[J]. IEEE Transactions on Visualization & Computer Graphics, 1999, 5(3):196-209.

[10]WAN M, DACHILLE F, KAUFMAN A. Distance-field based skeletons for virtual navigation[C]//Proceedings of the Conference on Visualization'01. IEEE Computer Society, 2001: 239-246.

[11]WANG S, WU J, WEI M, et al. Robust curve skeleton extraction for vascular structures[J]. Graphical Models, 2012, 74(4): 109-120.

[12]潘鵬, 賀三維, 吳艷蘭, 等. 曲邊多邊形中軸提取的新方法[J]. 測繪學(xué)報, 2012, 41(2): 278-283.

PAN Peng, HE Sanwei, WU Yanlan, et al. A new method for extracting curved-polygon medial axis[J]. Acta Geodaetica et Cartographica Sinica, 2012, 41(2): 278-283.

[13]CHEN W, SUI L, XU Z, et al. Improved Zhang-Suen thinning algorithm in binary line drawing applications[C]// International Conference on Systems and Informatic. Yantai:IEEE, 2012: 1947-1950.

[14]CHOI W P, LAM K M, SIU W C. Extraction of the Euclidean skeleton based on a connectivity criterion[J]. Pattern Recognition, 2003, 36(3): 721-729.

[15]MICHEL C, GILLES B. New Characterizations of Simple Points in 2D, 3D, and 4D Discrete Spaces[J]. IEEE Transactions on Pattern Analysis & Machine Intelligence, 2008, 31(4):637-648.

[16]GENG H, YANG J, TAN W, et al. Fast 3D skeleton extraction of airways and applications to virtual bronchoscopy[C]// The 26th Chinese Control and Decision Conference (2014 CCDC). Changsha:IEEE, 2014: 3879-3884.

[17]LIU L, CHAMBERS E W, LETSCHER D, et al. Extended grassfire transform on medial axes of 2D shapes[J]. Computer-Aided Design, 2011, 43(11): 1496-1505.

[18]GOLDBERGERA L, AMARAL L A N, GLASS L, et al. Physiobank, physiotoolkit, and physionet components of a new research resource for complex physiologic signals[J]. Circulation, 2000, 101(23): e215-e220.

[19]PALAGYI K. A 3D fully parallel surface-thinning algorithm[J]. Theoretical Computer Science, 2008, 406(1-2):119-135.

[20]楊威, 郭科, 魏義坤. 一種有效的基于八鄰域查表的指紋圖像細化算法[J]. 四川理工學(xué)院學(xué)報(自然科學(xué)版), 2008, 21(2): 61-63.

YANG Wei, GUO Ke, WEI Yikun. An effective thinning algorithm for finger print based on 8-neighbors LUT[J]. Journal of Sichuan University of Science and Engineering(Natural Science Edition), 2008, 21(2): 61-63.

[21]JIANG G, GU L. An automatic and fast centerline extraction algorithm for virtual colonoscopy[C]// 27th Annual International Conference of the Engineering in Medicine and Biology Society. Shanghai:IEEE, 2005: 5149-5152.

[22]陳磊, 王勝軍, 鄭全錄, 等. 基于 CT 圖像的三維拓撲細化算法及其在心臟 CAD 中的應(yīng)用[J]. 計算機應(yīng)用, 2007, 27(B06): 406-410.

CHEN Lei, WANG Shengjun, ZHENG Quanlu, et al. A 3D topological thinning algorithm based on CT image and its application in cardiac CAD, Computer Application, 2007, 27(B06): 406-410.

(編輯 王小唯， 苗秀芝)

A fast LUT-based airway skeleton extraction algorithm for virtual bronchoscopy

LIU Mingwei, GU Lixu

(School of Biomedical Engineering, Shanghai Jiao Tong University, Shanghai 200240, China)

In order to improve the efficiency of lung airway skeleton extraction, this paper introduces a new look-up-table (LUT) based thinning algorithm. This new approach consists of three major steps: the analysis and creation of the thinning model and the establishment of LUT based on above thinning model result. Thinning process by index-searching uses the LUT. Branch cutting based on the result of thinning obtains the final result. The usage of LUT index-searching during the thinning step transforms simple point judgments into LUT index search and has significantly improved the performance of the whole algorithm. Experimental results demonstrate that the new skeleton algorithm is 22.95 times faster than the existing thinning algorithm.

computer-aided diagnosis; lung airway; virtual bronchoscopy; skeleton; LUT

10.11918/j.issn.0367-6234.201601026

2016-01-01

國家自然科學(xué)基金(61271318)

劉明威(1990—)，男，碩士; 顧力栩(1964—)，男，教授，博士生導(dǎo)師

顧力翔，gulixu@sjtu.edu.cn.

文獻標(biāo)志碼: A 文章編號: 0367-6234(2017)05-0134-07