一種快速的視網(wǎng)膜光學(xué)相干圖像分割方法

張?zhí)鞓?黎日昌 羅文

摘要：在醫(yī)學(xué)領(lǐng)域，黃斑厚度可以用來(lái)量化糖尿病黃斑水腫和年齡相關(guān)性黃斑變性等疾病，臨床上通常使用光學(xué)相干斷層掃描的影像技術(shù)來(lái)獲取黃斑圖像。但現(xiàn)有的黃斑圖像分割方法運(yùn)算速度較慢，阻礙了其臨床使用。本文提出一種新的基于多分辨率及水平集的黃斑圖像分割方法，首先使用高斯濾波對(duì)原始圖像按行進(jìn)行濾波，再運(yùn)用多分辨率方法獲取圖像初始局部輪廓，最后使用水平集方法可以快速獲取黃斑圖像的中間輪廓，得到最終的圖像分割結(jié)果。通過(guò)在311幅黃斑圖像的仿真實(shí)驗(yàn)對(duì)比，本文方法在邊緣檢測(cè)結(jié)果和運(yùn)算速度上比傳統(tǒng)方法有很大改進(jìn)。

關(guān)鍵詞：黃斑圖像；光學(xué)相干成像；水平集；圖像分割

中圖分類(lèi)號(hào)：R770.43文獻(xiàn)標(biāo)識(shí)碼：A

1引言

光學(xué)相干斷層掃描（OCT）是一種可以通過(guò)光學(xué)信號(hào)采集和處理來(lái)對(duì)組織內(nèi)部結(jié)構(gòu)進(jìn)行無(wú)創(chuàng)傷成像的診斷設(shè)備。與其他高分辨率成像方法（如超聲波、X射線和磁共振成像）相比，光學(xué)相干斷層掃描方法具有更高的分辨率，可以顯示出更加清晰的組織內(nèi)部結(jié)構(gòu)圖像。目前，OCT這種光學(xué)成像技術(shù)已成為一種重要的可用于臨床上視網(wǎng)膜成像的實(shí)用技術(shù)。例如，臨床上通常使用光學(xué)相干斷層掃描的影像技術(shù)來(lái)獲取黃斑圖像。早期的OCT視網(wǎng)膜圖像分割方法主要是基于灰度閾值和灰度變化[1]，這些方法對(duì)噪聲敏感并且耗時(shí)。Koozekanani[2]等人提出了一種馬爾科夫隨機(jī)場(chǎng)（Markovrandomfield，MRF）的方法來(lái)提取視網(wǎng)膜內(nèi)部和外部邊緣，這種自回歸模型的魯棒性?xún)?yōu)于那些基于灰度閾值的方法。但是它需要可靠的初始種子點(diǎn)才能完成對(duì)病理視網(wǎng)膜的分割。Mujat[3]等人采用了形變樣條的方法來(lái)分割視網(wǎng)膜視神經(jīng)層，它需要把形變樣條放于初始輪廓的附近，通常也比較耗時(shí)。Chiu[4]指出文獻(xiàn)中報(bào)導(dǎo)的用于二維和三維OCT圖像分割的方法大多數(shù)速度較慢，這導(dǎo)致他們?cè)谂R床上實(shí)用性不足。

為了克服現(xiàn)有方法在速度上的不足，本文提出一種新的基于多分辨率及水平集的黃斑圖像分割方法。通過(guò)在311幅黃斑圖像的仿真實(shí)驗(yàn)對(duì)比，本文方法在邊緣檢測(cè)結(jié)果和運(yùn)算速度上比傳統(tǒng)方法有很大改進(jìn)。

2數(shù)據(jù)來(lái)源

通過(guò)光學(xué)相干斷層掃描成像方法，我們從34個(gè)病人中采集了總共311幅黃斑圖像，原始圖像的分辨率為2000（深度方向）*2048（寬度方向）。

3基于多分辨率及水平集的黃斑圖像分割方法

為了輔助醫(yī)學(xué)的黃斑厚度測(cè)量，需要得到黃斑圖像的清晰輪廓，本文設(shè)計(jì)一種新的基于多分辨率及水平集的黃斑圖像分割方法，首先使用一維高斯濾波對(duì)原始圖像按行進(jìn)行濾波，再運(yùn)用多分辨率方法獲取圖像初始局部輪廓，最后使用水平集方法可以快速獲取黃斑圖像的中間輪廓，得到最終的圖像分割結(jié)果。本方法的處理過(guò)程如圖1所示。

3.1高斯濾波

高斯濾波是一種線性平滑濾波，適用于消除高斯噪聲，廣泛應(yīng)用于圖像處理的減噪過(guò)程。

最常用的徑向基函數(shù)是高斯核函數(shù)，其函數(shù)形式為：

k（‖x-xc‖）=e-‖x-xc‖222σ2（1）

其中，xc為核函數(shù)中心，σ為核函數(shù)的寬度參數(shù)，決定了函數(shù)的徑向作用范圍。

本文使用高斯濾波對(duì)源圖像的每行進(jìn)行濾波。需要說(shuō)明我們只按行濾波，不在列方向上使用高斯濾波的原因是這樣將會(huì)導(dǎo)致邊緣位置的移動(dòng)。此外，我們使用一維高斯濾波，這是因?yàn)槎S高斯濾波速度較慢，而且可以分解為兩個(gè)一維高斯濾波。

3.2多分辨率分析

多分辨率分析是一種局部化時(shí)頻分析算法，可以用時(shí)域和頻域的聯(lián)合來(lái)表示信號(hào)，是分析非平穩(wěn)性信號(hào)的有力工具。它通過(guò)基函數(shù)的伸縮、平移等運(yùn)算多信號(hào)進(jìn)行多尺度細(xì)化分析，能有效從信號(hào)中提取信息，是一種靈活、快速、有效的高維信號(hào)處理算法。

本文使用多分辨率方法獲取黃斑圖像的初始輪廓，具體實(shí)現(xiàn)步驟包括以下步驟。首先，采用一個(gè)9×9的均值濾波器對(duì)高斯濾波后的黃斑圖像濾波，然后使用多分辨率的尺度H1W0，將圖像在豎直方向壓縮至源圖像的一半。為了提高分辨率，我們重復(fù)執(zhí)行若干次尺度變換，包括后續(xù)的尺度H2W0，尺度H3W1以及最后的尺度H4W2。

3.3水平集算法

水平集算法是Sethian和Osher[5]于1988年提出，用于把低維度的一些計(jì)算上升到更高一維，把N維的描述看成是N+1維的一個(gè)水平。

通常，水平集方法將平面閉合曲線隱含表達(dá)為連續(xù)函數(shù)曲面φ（x，y，t）的一個(gè)具有相同函數(shù)值的同值曲線。例如，可以將目標(biāo)曲線隱含表示在零水平集φ（x，y，t）=0中。即t時(shí)刻對(duì)應(yīng)于零水平集為

C（p，t）={（x，y）|φ（x，y，t）=0}（2）

假設(shè)用于演化的平面閉合曲線為C（p，t）=（x（p，t），y（p，t）），其中p為任意的參數(shù)化變量，t為時(shí)間。設(shè)曲線的內(nèi)向法向量為N，曲率為k，則曲線沿著其法向量的方向演化可以用下面的偏微分方程表示

Ct=V（k）N（3）

由于φ（C（t），t）=0，對(duì)t進(jìn)行全微分，整理得到水平集進(jìn)行曲線演化的方程為

SymbolQC@ φ（4）

從上述分析可知，水平集方法實(shí)現(xiàn)主動(dòng)輪廓線模型有幾個(gè)優(yōu)點(diǎn)：首先，演化曲線可以隨著φ的演化而改變拓?fù)浣Y(jié)構(gòu)，可以分裂、合并形成尖角等；其實(shí)，由于φ在演化過(guò)程中保持為一個(gè)完整的函數(shù)，容易實(shí)現(xiàn)近似數(shù)值計(jì)算；第三，水平集方法可以擴(kuò)展到高維曲面的演化，簡(jiǎn)化三維分割的復(fù)雜性。因此，我們可以采用快速水平集方法獲取中間輪廓。

具體實(shí)現(xiàn)過(guò)程中，我們采用了水平集方法來(lái)獲得黃斑層的邊緣輪廓，使用NumfordShah[6]水平集模型，它的基本公式如下：

其中，u和v是非負(fù)的，Ω表示圖像中的區(qū)域，曲線C表示曲線的邊緣，u表示初始圖像u0的片段。為了提高演化速度，我們采用了窄帶方法。

對(duì)于OCT黃斑圖像進(jìn)行水平集分割計(jì)算中間輪廓后，得到的黃斑層邊緣的中間分割結(jié)果如圖2所示。

圖2黃斑層水平集方法分割結(jié)果

3.4采用移動(dòng)多項(xiàng)式回歸的邊緣光順

從圖2可知，采用水平集方法獲取的邊緣是不光滑的。我們可以采用移動(dòng)多項(xiàng)式回歸的方法來(lái)對(duì)邊緣曲線進(jìn)行光順化，使其曲線光滑，去除噪聲。

多項(xiàng)式回歸方法

在本文中，我們采用移動(dòng)移動(dòng)多項(xiàng)式回歸的方法對(duì)邊緣曲線進(jìn)行光順。也就是以當(dāng)前點(diǎn)為中心，前后各取16個(gè)點(diǎn)，總共33個(gè)數(shù)據(jù)進(jìn)行多項(xiàng)式回歸，多項(xiàng)式的階數(shù)為5階。采用該方法后的結(jié)果如圖3.

4結(jié)論

針對(duì)現(xiàn)有黃斑圖像分割算法的邊緣提取存在著效率和速度的不足，本文提出一種新的基于多分辨率及水平集的黃斑圖像分割方法。首先使用一維高斯濾波對(duì)原始圖像按行進(jìn)行濾波，再運(yùn)用多分辨率方法獲取圖像初始局部輪廓，最后使用水平集方法可以快速獲取黃斑圖像的中間輪廓，得到最終的圖像分割結(jié)果。后續(xù)工作將在以下幾個(gè)方向上改進(jìn)：需要進(jìn)行參數(shù)和模型優(yōu)化，找到更有效的水平集或者主動(dòng)輪廓線算法來(lái)獲得邊緣輪廓，以進(jìn)一步提高圖像分割精度，其次還需要考慮如何分割不同病變圖像。

參考文獻(xiàn)

[1]CABRERAD.Delineatingfluidfilledregionboundariesinopticalcoherencetomographyimagesoftheretina[J].IEEETrans.Med.Image，2005，24，929-945.

[2]KOOZEKANANID，BOYER，KL，ROBERTSC.RetinalThicknessMeasurementsinOpticalCoherenceTomographyUsingaMarkovBoundaryModel[J].IEEETransactionsonMedicalImaging，2001，20（9），900-916

[3]MUJATM，CHANR，CENSEB，PARKB，JOOC，AKKINT，CHENT，deBOERJ.Retinalnervefiberlayerthicknessmapdeterminedfromopticalcoherencetomographyimages[J].OpticsExpress，2005，13（23），9480-9491

[4]CHIUS，LIC，NICHOLASP，TOTHC，IZATTJ，andSinaFarsiu.AutosegmentationofsevenretinallayersinSDOCTimagescongruentwithexpertmanualsegmentation[J].OpticsExpress，2010，18（18），19413-28.

[5]OSHERS，SETHIANJA.FrontsPropagatingwithCurvatureDependentSpeed：AlgorithmsBasedonHamiltonJacobiFormulations[J].JournalofComputationalPhysics，1988，79，12-49.

[6]MUMFORDD，SHAHJ.OptimalApproximationsbyPiecewiseSmoothFunctionsandAssociatedVariationalProblems[J].CommunicationsonPureandAppliedMathematics，1989，XLII（5），577-685.