郭志濤,周 峰,趙琳琳,袁金麗,盧成鋼
河北工業(yè)大學(xué) 電子信息工程學(xué)院,天津 300401
計(jì)算機(jī)斷層掃描(computed tomography,CT)技術(shù)能夠無侵入地對(duì)人體特定部位進(jìn)行快速精準(zhǔn)掃描而被廣泛應(yīng)用于臨床醫(yī)學(xué)放射診斷領(lǐng)域,然而CT掃描過程中產(chǎn)生的電離輻射會(huì)對(duì)人體造成傷害并增加患癌風(fēng)險(xiǎn)。目前,在實(shí)際臨床中一次CT掃描的輻射劑量在2~15 mSv之間,如一次胸部CT掃描的輻射劑量約為8 mSv,一次腹部CT掃描的輻射劑量約為10 mSv,而醫(yī)學(xué)上建議每年接受的輻射劑量不超過20 mSv[1],由于CT輻射具有累積效應(yīng),對(duì)于需要連續(xù)接受CT掃描的患者可能會(huì)造成嚴(yán)重的不良反應(yīng)。
對(duì)于如何改善、提高低劑量CT圖像(low-dose CT,LDCT)質(zhì)量的問題,現(xiàn)有的方法主要有三大類:投影域?yàn)V波算法,迭代重建算法和后處理方法[2]。其中全變分(TV)去噪是目前在這一領(lǐng)域得到成功應(yīng)用的最先進(jìn)的去噪方法之一。Zamyatin等[3]提出了一種自適應(yīng)的多尺度TV濾波算法,在保留圖像邊緣方面獲得了不錯(cuò)的效果。非局部(NLM)算法可以很好地利用圖像中自有的相似性和冗余性,Ma等[4]利用正常劑量掃描作為先驗(yàn)信息,將NLM方法應(yīng)用于常規(guī)劑量掃描中的信息冗余,在去除噪聲和抑制偽影方面效果顯著。為了提高圖像質(zhì)量,Chen等[5]采用了一種改進(jìn)的基于上下文的三維塊匹配算法(block matching 3D,BM3D),實(shí)現(xiàn)了較高質(zhì)量的去噪效果。字典學(xué)習(xí)和稀疏表示這一類方法也被廣泛應(yīng)用于LDCT圖像恢復(fù)。呂牡丹等[6]基于稀疏表示理論,引入了一種快速字典學(xué)習(xí)方法來處理低劑量CT圖像,有效地抑制了LDCT圖像中的斑點(diǎn)噪聲和條紋偽影。
近年來,隨著神經(jīng)網(wǎng)絡(luò)和深度學(xué)習(xí)在圖像處理領(lǐng)域的快速發(fā)展,Chen等[7]首次將卷積神經(jīng)網(wǎng)絡(luò)應(yīng)用于低劑量CT的后處理中,采用深度卷積神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu),將低劑量CT圖像逐塊映射到相應(yīng)的正常劑量CT(normaldose CT,NDCT)圖像,與傳統(tǒng)方法相比,在視覺效果與評(píng)價(jià)指標(biāo)上都顯示出一定的優(yōu)勢。Chen等[8]提出的具有殘差編解碼結(jié)構(gòu)的RED-CNN網(wǎng)絡(luò)在去噪效果上取得了優(yōu)異的結(jié)果,成為了低劑量CT去噪領(lǐng)域的經(jīng)典網(wǎng)絡(luò)結(jié)構(gòu)。Wu等[9]提出的級(jí)聯(lián)CNN去噪算法相較常規(guī)CNN算法而言,獲得了更高質(zhì)量的LDCT去噪圖像。Wolterink等[10]首先將生成對(duì)抗網(wǎng)絡(luò)(generative adversarial network,GAN)應(yīng)用于LDCT圖像去噪。Yang等[11]在GAN網(wǎng)絡(luò)中采用Wasserstein距離,形成了WGAN網(wǎng)絡(luò)來進(jìn)行低劑量CT圖像的去噪。Ge等[12]提出了ResvoxNet結(jié)構(gòu),將3D多尺度殘差結(jié)構(gòu)引入GAN網(wǎng)絡(luò)中,并采用成對(duì)的鑒別器來引導(dǎo)生成器復(fù)原出高質(zhì)量的LDCT圖像,取得了良好的效果。Yin等[13]基于漸近式三維殘差卷積網(wǎng)絡(luò)提出了DP-ResNet用于LDCT去噪。Chung等[14]提出了一種在空間頻率域和圖像域中運(yùn)行的雙域級(jí)聯(lián)U網(wǎng)絡(luò),以增強(qiáng)低劑量CT圖像的圖片質(zhì)量。
上述算法在去噪方面都取得了較好的效果,但是存在邊緣保護(hù)不足、圖像細(xì)節(jié)過于平滑的問題[15]。針對(duì)上述問題,本文提出了一種基于邊緣保護(hù)和多階段網(wǎng)絡(luò)的去噪模型。利用可訓(xùn)練的Sobel卷積,來提取輸入圖像的邊緣信息,實(shí)現(xiàn)模型對(duì)于邊緣細(xì)節(jié)的保護(hù)。同時(shí)采用多階段的去噪網(wǎng)絡(luò),逐級(jí)提取圖像特征,再通過特征融合整合各階段特征,從而得到更好輸出結(jié)果。
LDCT圖像具有較多噪聲,不利于臨床診斷。本文在圖像域通過深度學(xué)習(xí)方法去除LDCT圖像中的噪聲以恢復(fù)出高質(zhì)量的CT醫(yī)療圖像,LDCT圖像的恢復(fù)問題可以通過建立以下數(shù)學(xué)模型進(jìn)行解決。假設(shè)X∈Rm×n為低劑量CT圖像,Y∈Rm×n為對(duì)應(yīng)的正常劑量CT圖像,它們之間的關(guān)系可以表示為:
其中σ表示去噪過程。由LDCT圖像X恢復(fù)到NDCT圖像Y,即LDCT圖像的恢復(fù)過程可看作X=σ(Y)的逆映射,圖像重建的目標(biāo)是尋找將LDCT圖像X映射到NDCT圖像Y的函數(shù)F,使得如下所示的代價(jià)函數(shù)最小化:
其中F為σ-1的最佳近似值,在深度學(xué)習(xí)里面用可學(xué)習(xí)的卷積神經(jīng)網(wǎng)絡(luò)來表示F。
在邊緣保護(hù)模塊中,本文引入了可以提取邊緣特征的Sobel卷積,同時(shí)在Sobel算子中引入可訓(xùn)練的超參數(shù)α,使得算子可以在訓(xùn)練過程中不斷調(diào)整,從而提取合適的邊緣特征,如圖1(a)所示。此外,根據(jù)Sobel濾波的原理,將0°、45°、90°、135°四個(gè)方向上的運(yùn)算算子定義為一組。
圖1 Sobel算子示意圖Fig.1 Schematic image of Sobel operator
Sobel算子的運(yùn)算過程如圖1(b)所示,沿著兩個(gè)方向進(jìn)行運(yùn)算。在Sobel卷積的選擇過程中可以選擇多組Sobel算子來進(jìn)行運(yùn)算,在運(yùn)算結(jié)束之后,得到所提取的邊緣信息。
圖2的A為邊緣保護(hù)模塊(edge protect model,EPM)的實(shí)現(xiàn)過程,首先通過Sobel卷積,提取邊緣特征,再與輸入圖像在通道維度上進(jìn)行拼接。最后通過1×1卷積和Relu激活層得到最終的邊緣特征輸出。將邊緣信息與輸入圖像疊加起來,豐富了模型的輸入信息,加強(qiáng)了邊緣信息在網(wǎng)絡(luò)結(jié)構(gòu)中的作用。
如圖2所示,本文所提出的圖像重建算法框架包括三個(gè)階段來逐步恢復(fù)圖像。并引入邊緣保護(hù)模塊來實(shí)現(xiàn)對(duì)重建圖像邊緣細(xì)節(jié)的保護(hù)。
前兩個(gè)階段是基于編碼器-解碼器網(wǎng)絡(luò),采用較大的感受野,能夠?qū)W習(xí)廣泛的上下文信息[16]。為了使最終的重建圖像能夠保留更多圖像細(xì)節(jié),第三個(gè)階段采用能夠直接處理原始圖像的原始分辨率模塊(original resolution block,ORB),這一階段網(wǎng)絡(luò)主要由邊緣保護(hù)模塊、通道注意力模塊(channel attention block,CAB)以及卷積組成,不含有對(duì)圖像分辨率的操作,從而保留更多的圖像邊緣細(xì)節(jié)。
為了使不同階段的網(wǎng)絡(luò)可提取不同尺度的特征并且實(shí)現(xiàn)不同階段之間的特征融合,每個(gè)階段采用不同尺度的輸入,第一階段將輸入圖像分為四塊,即采用1/4的輸入尺度,第二階段采用1/2尺度作為輸入,最后一個(gè)階段采用原始圖像作為輸入。
(1)原始分辨率模塊:原始分辨率模塊實(shí)現(xiàn)過程如圖2的B所示,由通道注意力機(jī)制與1×1卷積構(gòu)成。該模塊直接作用于原始分辨率圖像,生成具有高分辨率與豐富圖像信息的特征,避免圖像信息的損失。
圖2 整體網(wǎng)絡(luò)結(jié)構(gòu)模型Fig.2 Overall network structure model
(2)通道注意力機(jī)制:本文所提方法通過通道注意力機(jī)制對(duì)提取的多通道特征進(jìn)行整合,如圖2的C所示。通過最大平均值池化生成特征圖,賦予標(biāo)記不同通道重要程度的權(quán)值。
(3)編解碼結(jié)構(gòu):該結(jié)構(gòu)基于一個(gè)三層標(biāo)準(zhǔn)的U-Net網(wǎng)絡(luò),如圖2中虛線方框區(qū)域所示,通過編碼器逐步地提取圖像的高維特征,再通過解碼器反向逐步恢復(fù)到輸入特征維度,使網(wǎng)絡(luò)可以提取到更多的上下文信息。
(4)特征融合過程:為了整合不同階段之間的特征,本文引入了跨階段特征融合過程。圖2的D為第一、二階段之間的特征融合過程,分為兩個(gè)部分:第一部分為由SAM提供的多通道特征;第二部分為上下兩個(gè)階段編碼器解碼器之間特征整合,通過1×1卷積來對(duì)特征進(jìn)行精細(xì)化,使網(wǎng)絡(luò)信息不容易丟失。圖2的E為二、三階段之間的特征融合過程。相比于一、二階段,加入了EPM模塊來保護(hù)邊緣細(xì)節(jié),加入特征的位置為ORB的輸入部分。
(5)監(jiān)督注意力模塊(supervise attention model,SAM)如圖3所示,在前兩個(gè)階段的最后輸出部分引入該模塊。首先,該模塊可以產(chǎn)生相應(yīng)階段的恢復(fù)圖像與正常劑量圖像的損失,為每個(gè)階段的圖像恢復(fù)提供監(jiān)督。其次,通過Sigmoid來生成注意力圖,對(duì)不同重要程度的特征進(jìn)行標(biāo)注,區(qū)分該階段不同信息量的特征,將有用的特征傳遞到下一階段進(jìn)行特征融合過程。
圖3 監(jiān)督注意力模塊(SAM)Fig.3 Supervise attention module(SAM)
SAM的輸入Fin表示由解碼器輸出的具有H×W×C維度的特征,本文先通過一個(gè)1×1卷積來生成圖像Rs,同時(shí)將維度為H×W×1的低劑量圖像與圖像Rs相加,生成該階段最終的預(yù)測圖像。再將預(yù)測圖像通過卷積和Sigmoid函數(shù)生成具有C通道的像素注意力掩模,使用這些掩模與經(jīng)過1×1卷積精細(xì)化之后的輸入特征Fin相乘,最后再與輸入特征相加產(chǎn)生注意力引導(dǎo)的特征。監(jiān)督注意力模塊的輸出Fout被傳遞到一下階段使用。
由于本文采用的為多階段網(wǎng)絡(luò)結(jié)構(gòu),將每個(gè)階段原始輸入圖像設(shè)為Xi,F(xiàn)為公式(2)中的映射函數(shù),Y為正常劑量圖像。α為邊緣保護(hù)模塊中引入的可訓(xùn)練參數(shù)。引入Charbonnier loss來優(yōu)化本文所提出的端到端網(wǎng)絡(luò)結(jié)構(gòu)。Charbonnier loss公式如下:
與MSE損失函數(shù)形式類似,其中ε表示一項(xiàng)極小項(xiàng),通常通過經(jīng)驗(yàn)來設(shè)定,在本文的實(shí)驗(yàn)中設(shè)置為10-3。Y表示NDCT作為標(biāo)簽來計(jì)算損失函數(shù)的值。
對(duì)于本文所提的多階段方法,本文綜合考慮各個(gè)階段的損失,用S來代表不同的階段,(本文所提結(jié)構(gòu)具有三個(gè)階段,故S為1,2,3)最終的損失函數(shù)可以表示為:
本章使用由Mayo Clinics授權(quán)的“2016 NIH-AAPMMayo Clinic低劑量CT大挑戰(zhàn)”的臨床數(shù)據(jù)集,用于網(wǎng)絡(luò)的訓(xùn)練和測試。數(shù)據(jù)集由投影數(shù)據(jù)和切片厚度分別為1 mm和3 mm重建圖像組成;使用3 mm厚的重建圖像,其中包含來自10位匿名患者的2 378張正常劑量CT圖像及其對(duì)應(yīng)的低劑量(1/4劑量)CT圖像,每張圖像的大小為512×512像素[17]。
實(shí)驗(yàn)中,選取來自9名患者的2 167對(duì)圖像作為訓(xùn)練集,另外一名患者的211對(duì)圖像作為測試集。本文的網(wǎng)絡(luò)為全卷積網(wǎng)絡(luò),因此能夠同時(shí)訓(xùn)練圖像塊和測試全尺寸圖像(512×512),從而有助于減少計(jì)算量,提高訓(xùn)練速度,并增加訓(xùn)練樣本數(shù),避免過擬合。
實(shí)驗(yàn)在64位Windows 10系統(tǒng)下進(jìn)行,采用PyTorch深度學(xué)習(xí)框架進(jìn)行網(wǎng)絡(luò)訓(xùn)練,并使用CUDA-Toolkit10.1對(duì)其加速;硬件配置為Intel Core i7-9700K CPU@3.2 GHz,內(nèi)存16 GB,顯卡為NVIDIA GeForce RTX 2080Ti。訓(xùn)練過程中使用Adam優(yōu)化算法對(duì)權(quán)值進(jìn)行迭代優(yōu)化。
2.3.1 主觀效果對(duì)比
為了綜合對(duì)比本文所提算法的性能,選取了如圖4所示的兩個(gè)樣本來進(jìn)行對(duì)比。如圖4(a)和圖4(b)為兩張腹部CT圖像,均具有較為復(fù)雜的圖像細(xì)節(jié)。
圖4 低劑量CT測試圖像Fig.4 Test images of low dose CT
本文選取的對(duì)比算法:BM3D算法、K-SVD算法、RED-CNN算法、WGAN-VGG算法。實(shí)驗(yàn)結(jié)果采用主觀效果和客觀指標(biāo)兩個(gè)方面進(jìn)行評(píng)價(jià)。視覺上以主觀效果突出圖像的結(jié)構(gòu)形態(tài),再利用客觀指標(biāo)評(píng)價(jià)算法的優(yōu)劣程度。
如圖5為圖4(a)的整體對(duì)比,圖6為圖5中ROI區(qū)域的放大對(duì)比。綜合圖5與圖6的對(duì)比,在視覺效果上看,所選用的對(duì)比算法都具備一定的去噪效果,如圖5(c)和圖5(e)分別為K-SVD和BM3D算法的結(jié)果圖,都去除了一定的噪聲,通過對(duì)比如圖5(f)、(d)和圖5(g)所示方法的結(jié)果圖可以看出,深度學(xué)習(xí)方法在視覺效果上相較于圖5(b)所示的LDCT圖像都有了一定提升,均顯示出了去噪的效果。通過對(duì)比可以看出,本文所提方法的恢復(fù)結(jié)果,如圖5(g)所示,在視覺效果上可以呈現(xiàn)出更多的圖像細(xì)節(jié)信息,圖像恢復(fù)結(jié)果更為清晰。在消除了噪聲和偽影的同時(shí)較好地保留了圖像的邊緣細(xì)節(jié),在成像質(zhì)量上與NDCT圖像最為接近。
圖5 測試樣本1整體對(duì)比圖Fig.5 Overall comparison of test sample 1
綜合圖6的ROI放大區(qū)域可以看出在箭頭所指部位邊緣細(xì)節(jié)過于模糊,損失了一些邊緣較為尖銳的細(xì)節(jié)。同時(shí)對(duì)比圖6的ROI放大圖來看,上述兩種去噪算法與NDCT圖像差距也較大。在整體視覺效果上看,本文所提方法在直觀的噪聲對(duì)比上具有更好的效果,與NDCT圖像也更為接近。通過對(duì)比圖6幾種去噪方法的ROI區(qū)域放大圖,可以實(shí)現(xiàn)對(duì)圖像邊緣的保護(hù),對(duì)比兩個(gè)箭頭所標(biāo)記出的部位,本文所提算法,圖6(g)在邊緣細(xì)節(jié)上更為清晰,與周圍背景的對(duì)比度也更為明顯。
圖6 測試樣本1 ROI區(qū)域放大對(duì)比圖Fig.6 Enlarged comparison of ROI area of test sample 1
如圖7為所選取測試樣本2的整體對(duì)比圖,圖8為圖7中ROI區(qū)域的放大對(duì)比圖。通過圖7的整體視覺效果對(duì)比可以看出,本文所提算法在不同結(jié)構(gòu)的CT圖像上也具有良好的去噪效果。對(duì)比圖8的ROI放大區(qū)域可以看出,本文所提算法不僅保留了較多的邊緣細(xì)節(jié),在直觀的去噪效果上也與NDCT圖像最為接近,整體畫面也較為干凈。
圖7 測試樣本2整體對(duì)比圖Fig.7 Overall comparison of test sample 2
通過圖5到圖8的對(duì)比結(jié)果可以看出,本文所提算法取得了最好的視覺效果,能夠保留與常規(guī)劑量圖像相近的細(xì)節(jié)信息,并除去了大部分噪聲,因此,從視覺效果看本文所提算法效果最佳。
圖8 測試樣本2 ROI區(qū)域放大對(duì)比圖Fig.8 Enlarged comparison of ROI area of test sample 2
2.3.2 客觀效果對(duì)比
為了客觀地評(píng)估本文所提算法的性能,本文采用峰值信噪比PSNR、空間結(jié)構(gòu)相似性SSIM作為低劑量CT圖像恢復(fù)效果的評(píng)價(jià)指標(biāo)來進(jìn)行圖像質(zhì)量評(píng)估。
如表1為測試樣本1和測試樣本2的客觀評(píng)價(jià)指標(biāo),其中最后一行為本文方法在測試樣本1、2中的客觀評(píng)價(jià)指標(biāo)結(jié)果??梢钥闯?,本文方法在兩張測試圖中均取得了最高的峰值信噪比和空間結(jié)構(gòu)相似性。這表明在這兩張測試圖中,本文所提算法在去噪效果方面和與NDCT圖像結(jié)構(gòu)相似性方面均取得了最好的結(jié)果。
表1 測試樣本1和測試樣本2客觀指標(biāo)Table 1 Objective indicators of test sample 1 and test sample 2
如表2為不同算法在測試集中的平均客觀指標(biāo)與預(yù)測時(shí)長對(duì)比,最優(yōu)結(jié)果用加黑字體標(biāo)出。可以看出,在客觀評(píng)價(jià)指標(biāo)PSNR和SSIM上,本文所提方法也取得了最好的結(jié)果,表明了本文算法在測試集中獲得最好的去噪效果同時(shí)也保留了與NDCT圖像最好的結(jié)構(gòu)相似性。
表2在平均預(yù)測每張圖片所用的時(shí)間的對(duì)比中,為了避免不同情況下預(yù)測時(shí)間的差異,本文對(duì)同一測試集分別進(jìn)行10次預(yù)測,并取預(yù)測平均值進(jìn)行了對(duì)比??梢钥闯?,深度學(xué)習(xí)方法的平均預(yù)測時(shí)間均要好于BM3D與K-SVD算法。本文所提算法的預(yù)測時(shí)長在優(yōu)于WGAN-VGG算法的情況下取得了更好的性能,生成的圖片更為精細(xì),保留了更多的圖像細(xì)節(jié)。同時(shí)本文所提算法與RED-CNN網(wǎng)絡(luò)的預(yù)測時(shí)長幾乎相同,且取得了更好的客觀評(píng)價(jià)指標(biāo)PSNR、SSIM與主觀視覺效果。
表2 測試集中不同算法客觀指標(biāo)與預(yù)測時(shí)長對(duì)比Table 2 Comparison of objective indexes and prediction duration of different algorithms in test set
將本文算法與最新的低劑量CT圖像去噪算法進(jìn)行了對(duì)比,如表3所示??梢钥闯鱿啾扔赗esvoxNet、DP-Resnet、FI-Wnet,本文所提方法在峰值信噪比和空間結(jié)構(gòu)相似性上均取得了最優(yōu)的評(píng)價(jià)指標(biāo)。
表3 不同算法的客觀指標(biāo)對(duì)比Table 3 Comparison of objective indexes of different algorithms
本文采用消融實(shí)驗(yàn),用來分析模型不同的組成部分所起到的作用,該部分實(shí)驗(yàn)與對(duì)比實(shí)驗(yàn)所選用數(shù)據(jù)集相同。
在消融實(shí)驗(yàn)中,主要通過對(duì)比邊緣保護(hù)模塊的引入與損失函數(shù)的選取來評(píng)估本文所提模型的性能。本文方法引入了邊緣保護(hù)模塊,由Sobel算子來實(shí)現(xiàn),對(duì)比了MSD-Net(multi stage denosing-Net)和所提算法。由表4前兩行與后兩行對(duì)比可以看出引入了邊緣保護(hù)模塊之后,PSNR、SSIM在指標(biāo)效果上均有所提升。本文同時(shí)也采用了疊加形式的Charbonnier loss(在表中表示為Closs)綜合考慮各階段的輸出對(duì)于最終去噪結(jié)果的影響。對(duì)比表中使用Closs和MSE loss的性能指標(biāo),可以看出本文所提方法取得了最優(yōu)的評(píng)價(jià)指標(biāo)。
表4 消融實(shí)驗(yàn)客觀指標(biāo)對(duì)比Table 4 Comparison of objective indexes of ablation experiments
對(duì)于低劑量CT圖像成像質(zhì)量低,現(xiàn)有后處理方法在重建過程中存在圖像過于平滑,邊緣細(xì)節(jié)不足的問題,本文提出了將傳統(tǒng)Sobel算子和近年來興起的深度學(xué)習(xí)相結(jié)合的方法。同時(shí)在綜合考慮各個(gè)階段對(duì)于最終恢復(fù)圖像的影響下引入了Charbonnier loss。本文方法充分利用了網(wǎng)絡(luò)空間中的上下文聯(lián)系,有效獲取了空間特征的全局信息,實(shí)現(xiàn)了網(wǎng)絡(luò)對(duì)關(guān)鍵特征位置的全面把握,提升了網(wǎng)絡(luò)輸出的信息的豐富度。在主觀視覺效果對(duì)比上,該網(wǎng)絡(luò)在保持圖像信息的同時(shí),能夠很好地降低圖像噪聲,保持器官的邊緣細(xì)節(jié)。對(duì)比客觀評(píng)價(jià)指標(biāo),本文方法在峰值信噪比(PSNR)和空間結(jié)構(gòu)相似度(SSIM)上優(yōu)于其他主流方法,同時(shí)在預(yù)測時(shí)長上也取得了優(yōu)秀的表現(xiàn)。綜上所述,本文所提算法在從低劑量CT圖像中恢復(fù)出正常劑量CT圖像方面具有一定的先進(jìn)性。