段繼忠, 錢青青

(昆明理工大學 信息工程與自動化學院,昆明 650504)

磁共振成像(Magnetic Resonance Imaging, MRI)利用人體組織中氫原子核的核磁共振現(xiàn)象,通過接收線圈采集k空間(即頻域)數(shù)據(jù),然后對采集到的k空間數(shù)據(jù)進行傅里葉逆變換,從而重建出人體器官圖像[1].MRI具有無電離輻射、多角度成像、對人體組織無損傷等諸多優(yōu)點,是臨床醫(yī)學和醫(yī)學科研中非常重要的檢測手段[2].

然而,MRI掃描時間較長,病人在掃描過程中易感到不適或發(fā)生移動導致圖像偽影,影響病人的就醫(yī)體驗和掃描圖像質(zhì)量.為加快MRI的掃描速度,通常對k空間數(shù)據(jù)進行欠采樣后再使用算法重建出圖像,但欠采樣易使圖像質(zhì)量下降,影響醫(yī)生對患者病情的診斷.在算法模型中加入不同正則項可改善圖像重建質(zhì)量,但會導致重建速度較慢從而出診斷結果較慢,最終影響病人的就醫(yī)滿意度,限制了算法的推廣與應用[3].因此,在保證MRI重建質(zhì)量的同時,研究MRI的快速重建方法是非常實際且有意義的.

并行成像技術[4]和壓縮感知(Compressed Sensing, CS)[5]是減少MRI掃描時間的兩種有效方法.其中,并行磁共振成像(Parallel MRI, PMRI)方法使用具有不同空間靈敏度信息的多通道接收線圈陣列同時采集磁共振信號,利用通道間的關聯(lián)信息即可重建出未采集的信號.而CS理論突破了奈奎斯特采樣定理,利用圖像的可壓縮性和稀疏性即可從很少的觀測值中重建出高質(zhì)量圖像.為進一步提高MRI的重建性能,通常將PMRI方法與CS理論結合應用[6].

在PMRI技術中,Lustig等[7]提出一種針對任意k空間軌跡的并行成像重建方法,即迭代自一致性并行成像重建(Iterative Self-Consistent Parallel Imaging Reconstruction, SPIRiT)模型.該方法基于k空間數(shù)據(jù)自身的一致性,包括校正一致性和數(shù)據(jù)采集一致性,利用相鄰k空間點之間的相關性恢復丟失的信息,已被廣泛應用于臨床實踐中[8].SPIRiT是一種基于自一致性的廣義PMRI重建框架,可以方便地與各種正則化方法結合,以實現(xiàn)更高性能的MRI重建.因此,許多研究者引入不同正則項對SPIRiT模型進行一系列的研究和改進[9-17].其中,段繼忠等[9]將聯(lián)合稀疏性正則項引入SPIRiT模型,提出一種高效的重建方法.之后,Duan等[12]又在SPIRiT模型中引入聯(lián)合全變分和聯(lián)合L1范數(shù)正則項約束,提出另一種基于SPIRiT的快速重建方法.這兩種方法都將復雜優(yōu)化問題進行簡化,然后利用算子分裂技術將其分解為易于計算或求解的子問題,最后使用快速迭代軟閾值算法(Fast Iterative Shrinkage-Thresholding Algorithm, FISTA)[18]進行加速,縮短圖像的重建時間.最近,Zhang等[15]將L1范數(shù)正則項加入SPIRiT模型中,提出一種基于快速投影迭代軟閾值算法(projected FISTA, pFISTA)[19]的快速重建方法(pFISTA-SPIRiT).該方法允許在并行MR圖像重建中使用不同的緊框架,是一個廣義模型,已被用于人腦成像、肝臟動態(tài)增強成像和精準的腦腫瘤血管通透性成像等[20].

但pFISTA-SPIRiT方法對校正一致項與數(shù)據(jù)一致項都進行梯度計算,導致收斂速度較慢.將重建模型中的復雜優(yōu)化問題進行簡化可以有效提升算法的收斂速度[9,12].因此,本文基于pFISTA-SPIRiT模型,提出一種基于平移不變離散小波變換(Shift-Invariant Discrete Wavelets Transform, SIDWT)和SPIRiT的快速并行磁共振成像重建方法(fast SIDWT-SPIRiT, fSIDWT-SPIRiT).與pFISTA-SPIRiT利用圖像域SPIRiT模型進行重建不同,該方法利用頻域SPIRiT模型進行重建.使用SIDWT和SPIRiT模型,將校正一致項和數(shù)據(jù)一致項合并為一項后再使用pFISTA技術進行求解.實驗結果表明,與其他兩種方法相比,該方法能夠有效提升算法的收斂速度,從而減少圖像的重建時間,并且重建質(zhì)量不變.

1 SPIRiT模型

SPIRiT是一種逐線圈、自動校準的并行MRI重建模型.設r為k空間數(shù)據(jù)點的位置,Rr為從k空間中選擇需要的數(shù)據(jù)點的算子;xj為第j個線圈上的全部k空間數(shù)據(jù),j=1,2,…,J,J為總線圈數(shù);Rrxj為以第j個線圈位置r為中心從所有線圈提取的k空間鄰域數(shù)據(jù)點;xi(r)為第i個線圈上位置r處的k空間數(shù)據(jù),則xi(r)的重建如下:

(1)

(2)

式中:A為ACS的位置集合.式(1)是一組耦合線性方程,故可以將整個耦合方程組寫成矩陣形式,則所有線圈的校正一致性如下:

x=Gx

(3)

式中:x為從所有線圈中獲取的全部k空間數(shù)據(jù);G為在適當位置包含gji的矩陣.

當然,除了考慮校正一致性,還應滿足數(shù)據(jù)采集一致性.數(shù)據(jù)采集一致性如下:

y=Dx

(4)

式中:y為從所有線圈中獲取的欠采樣數(shù)據(jù);D為從全k空間數(shù)據(jù)進行欠采樣的矩陣.

2 fSIDWT-SPIRiT

在重建模型中加入正則項可以有效提升MRI的重建質(zhì)量,因此在SPIRiT模型中引入L1范數(shù)正則項,則帶L1范數(shù)正則項的SPIRiT并行MRI重建模型可以表示為以下最優(yōu)化問題:

(5)

式中:Ψ為逐線圈小波變換,用于將圖像稀疏化,是一個緊框架.本文選用SIDWT作為實驗中的緊框架.F為逐線圈傅里葉變換,F-1為逐線圈傅里葉逆變換.

使用罰函數(shù)技術,可以得到式(5)的無約束版本:

(6)

(7)

針對式(7)中無約束模型的求解,令α=ΨF-1x,并且有Ψ*Ψ=I,則

(8)

(9)

式中:k為迭代次數(shù);L為f(αk)梯度的Lipschitz常數(shù).于是,易得式(9)的解為

(10)

(11)

式中:β為輸入矩陣;|β|為β的模.由于Ψ*Ψ=I故式(10)可重寫為

(12)

(13)

最后,利用FISTA[18]技術進行加速,則整個求解過程如下:

(14)

(15)

(16)

(17)

綜上所述,得到基于SIDWT和SPIRiT模型的快速并行磁共振成像重建方法fSIDWT-SPIRiT,實現(xiàn)步驟如下.

算法一: fSIDWT-SPIRiT

1: Setx0=0,z0=0,θ0=0,k=0

2: Repeat

3:k=k+1

7: 直到達到最大迭代次數(shù),否則返回步驟3

3 實驗仿真及分析

3.1 實驗設置

所有實驗均在計算機上進行,計算機的配置為:Intel(R) Core(TM) i5-7200U CPU@2.50 GHz處理器、16 GB內(nèi)存和Windows 10操作系統(tǒng)(64位),所有方法均使用MATLAB實現(xiàn).

在實驗中,將fSIDWT-SPIRiT的重建性能與pFISTA-SPIRiT和SIDWT-SPIRiT進行對比.后兩種方法都基于SIDWT和SPIRiT模型,直接使用pFISTA技術進行求解得到,區(qū)別在于pFISTA-SPIRiT在圖像域SPIRiT進行重建,而SIDWT-SPIRiT在k空間域SPIRiT進行重建.使用信噪比(Signal Noise Ratio, SNR)、結構相似性(Structural SIMilarity, SSIM)指標和高頻誤差范數(shù)(High-Frequency Error Norm, HFEN)3個評價指標來衡量圖像的重建質(zhì)量.3個評價指標的定義分別如下:

(18)

(19)

(20)

式中:VHFEN為HFEN值;filter(·)是一個拉普拉斯高斯濾波器,用于捕捉圖像邊緣,濾波核的大小為15×15像素,標準差為1.5像素.

VSNR和VSSIM的數(shù)值越高,VHFEN的數(shù)值越低,說明圖像的重建質(zhì)量越好.這3個評價指標均在感興趣區(qū)域內(nèi)進行計算,實驗中手動調(diào)整所有方法的參數(shù)使得VSNR值達到最優(yōu).

為驗證fSIDWT-SPIRiT的有效性,在不同人體器官數(shù)據(jù)集上對各種方法的重建性能進行比較,分別選取活體人類腦部切片圖GE_human_brain[21]、train_brain_AXT1POST_200_6001959[22-23]、心臟切片圖data_v1_k1[24]以及人類膝蓋切片圖train_knee_file1000005[22-23],并依次命名為數(shù)據(jù)集1,數(shù)據(jù)集2,數(shù)據(jù)集3和數(shù)據(jù)集4.其中,數(shù)據(jù)集1是尺寸為256×256像素的8通道腦部切片數(shù)據(jù);數(shù)據(jù)集2是使用20通道線圈獲取的腦部數(shù)據(jù)集的第1幀,然后使用線圈壓縮技術將其壓縮為8個虛擬線圈,尺寸為320像素×320像素;數(shù)據(jù)集3是使用28通道線圈獲取的心臟數(shù)據(jù)集,然后使用線圈壓縮技術將其壓縮為12個虛擬線圈,尺寸為192像素×192像素;數(shù)據(jù)集4是在15通道線圈獲取的膝蓋數(shù)據(jù)集的第20個切片,然后使用線圈壓縮技術將其壓縮為8個虛擬線圈,尺寸為320像素×320像素.

實驗中的所有數(shù)據(jù)集都使用具有不同加速因子R(不包括ACS)的二維泊松圓盤采樣模式進行欠采樣,并且所有方法都使用大小為24像素×24像素的校準區(qū)域和5像素×5像素的SPIRiT核.本文主要從重建質(zhì)量和重建速度兩個方面來對比各種方法的優(yōu)劣.在重建質(zhì)量方面,主要通過重建圖和誤差圖從主觀上比較各種方法的重建性能,并采用3個評價指標從客觀上比較各種方法的重建性能;在重建速度方面,主要通過時間來比較各種方法的重建性能.

3.2 不同方法的重建質(zhì)量比較

3.2.1不同重建方法的視覺比較 首先,從視覺上比較各種方法的重建性能.當加速因子R=5(即欠采樣率為20%)時,用3種方法對4個數(shù)據(jù)集進行重建.圖1～4分別給出4個數(shù)據(jù)集對應的全采樣圖像、二維泊松圓盤欠采樣掩膜、重建圖及誤差圖.

圖1 在5倍加速時不同方法對數(shù)據(jù)集1進行重建的視覺比較Fig.1 Comparison of visual reconstruction for dataset 1 by different methods at 5 times acceleration

圖2 在5倍加速時不同方法對數(shù)據(jù)集2進行重建的視覺比較Fig.2 Comparison of visual reconstruction for dataset 2 by different methods at 5 times acceleration

由圖1～4可以看出,對于4個不同人體器官的數(shù)據(jù)集,fSIDWT-SPIRiT與其余兩種方法相比,重建圖和誤差圖都無明顯差異.其中,fSIDWT-SPIRiT在不同數(shù)據(jù)集上的重建圖都有較清楚的紋理細節(jié),并保留了較好的邊緣輪廓信息,表明所提方法能夠?qū)崿F(xiàn)較好的重建,與其余方法的重建質(zhì)量相當.

3.2.2不同重建方法的評價指標比較 從3個評價指標上比較各種方法的重建性能.在加速因子R為3～7的情況下,即欠采樣率分別為33.3%、25%、20%、16.7%和14.3%時,用3種方法對4個不同數(shù)據(jù)集進行重建.表1～4分別給出4個不同數(shù)據(jù)集下不同方法對應的VSNR值、VSSIM值和VHFEN值.

表1 在3～7倍加速時不同方法對數(shù)據(jù)集1進行重建的數(shù)值比較Tab.1 Comparison of the values of the reconstruction for dataset 1 by different methods at 3 to 7 times acceleration

表2 在3～7倍加速時不同方法對數(shù)據(jù)集2進行重建的數(shù)值比較Tab.2 Comparison of the values of the reconstruction for dataset 2 by different methods at 3 to 7 times acceleration

表3 在3～7倍加速時不同方法對數(shù)據(jù)集3進行重建的數(shù)值比較Tab.3 Comparison of the values of the reconstruction for dataset 3 by different methods at 3 to 7 times acceleration

表4 在3～7倍加速時不同方法對數(shù)據(jù)集4進行重建的數(shù)值比較Tab.4 Comparison of the values of the reconstruction for dataset 4 by different methods at 3 to 7 times acceleration

從表1～4可以看出,對于數(shù)據(jù)集1和數(shù)據(jù)集3,與pFISTA-SPIRiT相比,fSIDWT-SPIRiT的重建圖像VSNR值和VSSIM值更高,VHFEN值更低,說明所提方法的重建質(zhì)量更好;而對于數(shù)據(jù)集2和數(shù)據(jù)集4,pFISTA-SPIRiT和fSIDWT-SPIRiT的重建圖像VSNR值、VSSIM值和VHFEN值相差并不大,說明兩種方法的重建質(zhì)量無顯著差異. 對于4個不同的數(shù)據(jù)集,fSIDWT-SPIRiT與SIDWT-SPIRiT相比,重建圖像的VSNR值、VSSIM值和VHFEN值都很接近,說明兩種方法的重建質(zhì)量相當.

3.3 不同方法的重建速度比較

最后,從速度上比較各種方法的重建性能.在加速因子R為3～7的情況下,用3種方法對4個數(shù)據(jù)集進行重建.圖5～6為R=3和R=5時3種方法對4個不同數(shù)據(jù)集進行重建的時間(t)比較.表5給出R為3～7時在4個不同數(shù)據(jù)集下fSIDWT-SPIRiT相對于其他方法重建時間提升的倍數(shù)(m).

表5 在3～7倍加速時不同方法對4個數(shù)據(jù)集進行重建的時間比較Tab.5 Comparison of reconstruction time for 4 datasets by different methods at 3 to 7 times acceleration

圖5 在3倍加速時3種方法對4個數(shù)據(jù)集的重建速度比較Fig.5 Comparison of the reconstruction speed of the three methods for the 4 datasets at 3 times acceleration

圖6 在5倍加速時3種方法對4個數(shù)據(jù)集的重建速度比較Fig.6 Comparison of the reconstruction speed of the three methods for the 4 datasets at 5 times acceleration

從圖5～6可以看出,對于4個不同的數(shù)據(jù)集,fSIDWT-SPIRiT與pFISTA-SPIRiT和SIDWT-SPIRiT兩種方法相比,無論加速因子R=3還是R=5,均能獲得比較接近的VSNR值.但fSIDWT-SPIRiT的收斂速度明顯快于其他兩種方法,表明所提方法可以實現(xiàn)更快速度的重建,并保證圖像的重建質(zhì)量,進一步驗證方法的有效性.

從表5可以看出,從3倍加速增加到7倍加速,對于4個不同的數(shù)據(jù)集,與pFISTA-SPIRiT相比,fSIDWT-SPIRiT的重建速度分別提升3.2倍、3.7倍、3.8倍和3.4倍,平均提升3.5倍;與SIDWT-SPIRiT相比,fSIDWT-SPIRiT的重建速度分別提升3.5倍、3.4倍、4.9倍和4倍,平均提升3.9倍.由此可以說明,所提方法能夠大大縮短圖像重建時間,更進一步驗證該方法的有效性.

3.4 方法的收斂性分析

探討VSNR值、VHFEN值以及相對誤差(Relative Error, RE)隨迭代次數(shù)k的變化,分析fSIDWT-SPIRiT的收斂性.RE的值計算如下:

(21)

如圖7所示,在3倍加速下,對于數(shù)據(jù)集1,在前10次迭代中VSNR值和VHFEN值變化很快,第10次迭代之后變化較為緩慢,當VRE<2×10-3即迭代次數(shù)大于16次時,VSNR值和VHFEN值基本不再變化.此時認為所提方法收斂,利用此特性提前終止算法可以減少計算量.對于其他3個數(shù)據(jù)集可以得出類似的結論.對于數(shù)據(jù)集2～4,當VRE分別小于5×10-3、3×10-3和4×10-3即迭代次數(shù)分別大于13次、23次和23次時,認為所提方法收斂,可以提前終止算法以減少計算量.

圖7 在3倍加速時對4個數(shù)據(jù)集進行重建的收斂性分析Fig.7 Convergence analysis of the reconstruction for 4 datasets at 3 times acceleration

4 結語

基于SIDWT和SPIRiT模型,利用pFISTA技術進行求解,提出一種新的快速并行磁共振成像重建方法——fSIDWT-SPIRiT.在4個不同活體數(shù)據(jù)集上的仿真實驗表明:與pFISTA-SPIRiT和SIDWT-SPIRiT兩種方法相比,所提方法能獲得與之相當?shù)闹亟ㄙ|(zhì)量,而且收斂速度明顯更快,平均提升3.5倍和3.9倍.因此,所提方法既保證圖像的重建質(zhì)量,還顯著減少圖像的重建時間.本文只選用緊框架SIDWT進行實驗,后續(xù)研究中將考慮使用其他緊框架進一步提升圖像的重建性能.