基于動態(tài)卷積的體素內不相干運動成像參數估計①

鄧開軒, 程欣宇, 王麗會

(貴州大學 計算機科學與技術學院 貴州省智能醫(yī)學影像分析與精準診斷重點實驗室, 貴陽 550025)

擴散加權磁共振(diffusion-weighted magnetic resonance imaging, DW-MRI)是一種以非侵入性的方式, 測量活體組織內水分子擴散運動特性的成像技術.在采集過程中, 通過改變擴散梯度強度與持續(xù)時間 (即擴散敏感因子,b值)來改變對水分子擴散運動的敏感性, 從而測定出衰減的磁共振信號. 通常利用單指數模型來衡量衰減的磁共振信號, 估計出表觀擴散系數(apparent diffusion coefficient, ADC), 進而間接地反映組織細胞的微觀結構[1].

體素內不相干運動成像(intravoxel incoherent motion imaging, IVIM)是一種常用的多b值DW-MRI技術,除了能夠提供人體組織間質中水分子運動的擴散信息之外, 還能夠提供源于血液微循環(huán)中的灌注信息[2].IVIM模型利用微循環(huán)在低b值會導致信號快速衰減這一特性, 通過建立關于b值與衰減信號的雙指數模型來分離擴散和灌注運動, 估計出灌注體積分數(fraction of prefusion,fp)、灌注或偽擴散系數(perfusion or pseudo-diffusion coefficient,Dp)和擴散系數(diffusion coefficient,Dt). 通過測定組織的擴散和灌注信息, IVIM已經在肝臟[3]、腎臟[4]、乳腺[5]、心臟[6]以及頭頸部[7]等多種灌注信息豐富的器官上用于臨床研究, 在疾病早期發(fā)現及術后評價中發(fā)揮重要的作用[8].

盡管IVIM成像技術在臨床實踐中取得了重大進展, 但從獲得的DW信號中快速準確地估計出高質量的灌注和擴散相關參數仍然是具有挑戰(zhàn)性的問題. 雖然在理論上更多的b值有利于提升IVIM參數估計的精度, 但是這無疑會使得圖像采集時間過長, 同時過多的低b值會導致圖像的信噪比(signal to noise ratio,SNR)降低. 因此解決這些缺點將有助于更加廣泛地使用IVIM. Zhang等提出了一種最小化誤差傳播因子的優(yōu)化b值采樣的方法, 顯著提高了通過非均勻優(yōu)化采樣獲得的腎臟IVIM參數準確度, 降低了至少20%–30%的估計誤差, 并顯著提高了惡性和良性病變的區(qū)別[9].Huang提出了一種基于核的后濾波圖像降噪方法, 與局部主成分分析法進行對比, 得到了均小于12%的偏差和變異系數的IVIM參數[10]. Lin等利用深度圖像先驗的思想(deep image prior, DIP)提出了一種無監(jiān)督的單圖像降噪網絡, 利用卷積神經網絡(convolutional neural network, CNN), 將多組在b值上平均IVIM數據作為輸入, 輸出與原始數據b值相同數量的降噪圖像,提高了DW數據的信噪比, 得到的IVIM參數圖干凈平滑[11]. Huang 等利用圖像在多尺度上層內與層間的多相似性, 通過隨機搜索策略將低分辨率圖像映射到高分辨率圖像空間中, 使得估計出的IVIM參數圖像有更低的平均殘差值[12]. 上述方法中, 分別從數據優(yōu)化采集策略、降低生理噪聲、圖像超分辨率等不同角度的數據預處理方法使IVIM參數質量得到了間接改善, 但是參數估計的準確性仍然依賴于擬合方法.

目前, IVIM最常用的擬合方法是非線性最小二乘法(nonlinear least squares method, NLLS), 即通過最小化模型的擬合信號與采集信號間的誤差來優(yōu)化參數,其中最常用的是LM (Levenberg Marquardt)優(yōu)化算法[13].這種方法受到了來自于采集信號信噪比嚴重的影響.為了提高NLLS擬合的穩(wěn)定性, 一般嘗試在擬合過程中根據具體估計的器官, 加入各個參數的約束項, 結果證明在參數約束范圍內進行擬合, 參數的評估質量得到了提升. 但是, 不少研究表明, 利用帶有約束的NLLS方法估計出來的IVIM參數常出現在約束邊緣值附近,造成了結果的偏離性. 根據組織內擴散和灌注信號衰減的特點, 通常采用分段式NLLS方法(segmented NLLS, Seg-NLLS)來提高參數估計的魯棒性[13–15]. Seg-NLLS首先考慮在高b值下, 組織灌注信號的貢獻為0,將雙指數模型簡化為單指數模型, 對擴散系數Dt進行估計. 接著利用估計好的擴散系數Dt, 推算出灌注分數fp和灌注系數Dp. 然而NLLS方法受噪聲影響嚴重,特別是灌注系數Dp, 在某些局部區(qū)域, 由于噪聲的影響將導致其估計結果無意義[16–20]. Hyo在12名健康志愿者的肝臟和胰腺數據進行了可重復性測試, 實驗結果表明上述基于最小二乘法得到的Dp的可重復性最差[16].為了解決上述問題, 陸續(xù)有學者采用概率的思想來評估IVIM參數, 即貝葉斯概率理論(Bayesian approach,Bayes). Bayes是一種以概率先驗替代最小二乘的IVIM參數求解方法, 它通過最大化給定觀測信號的IVIM參數的后驗概率來估計IVIM參數[21–23]. 最近一些工作比較了各IVIM估計方法的性能, 實驗結果表明Bayes方法始終優(yōu)于NLLS[18,20,23–25]. 然而, 雖然Bayes擁有較高的參數估計精確度, 但是該方法很容易受到組織異質性的影響, 并且估計精度強烈依賴于先驗假設, 擬合時間較長. 此外, 由于以參數為優(yōu)化目標, 其擬合信號與測量信號之間的誤差大小難以控制[26]. 為了解決傳統IVIM參數擬合方法的缺陷, 部分學者利用人工神經網絡(artificial neural network, ANN)對IVIM參數進行估計[26–29], ANN由輸出層、隱含層和輸出層構成, 本質是通過高維數據的低維非線性映射, 實驗結果表明ANN產生的IVIM參數圖具有更好的視覺效果和較高的一致性. 然而, ANN并沒有考慮到到信號之間的上下文信息, 同時ANN也容易過擬合. 針對于上述IVIM參數方法的不足, 本文提出一維動態(tài)卷積網絡(DCNN),考慮信號之間的上下文信息和b值對IVIM參數的貢獻程度, 實現IVIM參數的準確估計. 仿真和真實采集兩種數據集的實驗表明, DCNN有更好的IVIM參數估計性能和更好的視覺質量.

1 相關技術與對比方法

1.1 體素內不相干運動模型

為了評估這些擬合方法的IVIM參數估計, 我們使用帶噪聲的雙指數模型來測定水分子在人體大腦組織中的擴散信息和灌注信息, 該模型將DW-MRI衰減信號描述成關于b值即擴散敏感度(diffusion sensitivity,b-value)的函數. 在給定的b值下, IVIM公式如式(1):

其中,n表示b值個數,S(bn)表示為第n個b值下的擴散加權信號,S(b0)表示為b= 0 s/mm2時的初始信號強度.Dt為擴散系數,Dp為灌注系數,fp為灌注體積分數,εn表示信號含有的噪聲.

1.2 非線性最小二乘法

NLLS通過最小化模型的擬合信號與采集信號間的誤差平方和來優(yōu)化參數, 如式(2)所示:

其中,Sint代表在高b值下數據的截距. 灌注體積分數fp再由式(4)所示計算出:

接下來, 在得到體素點的Dt和fp后, 再通過一次最小二乘擬合出Dp. 由于Seg-NLLS會導致偏置誤差,在本研究中將使用Full-NLLS作為基礎參數擬合算法,當體素點對應的IVIM參數超過擬合邊界時, 再使用Seg-NLLS對該體素點進行估計.

1.3 貝葉斯概率法

Bayes是通過最大化IVIM參數的后驗概率密度函數來估計, 并使用馬爾科夫鏈蒙特卡羅方法(Markov chain Monte Carlo, MCMC)實現[30]. 多參數的高斯似然函數如式(3)所示:

其中,Dt和Dp采用對數正態(tài)分布先驗,fp采用beta分布先驗,代表高斯收縮的方差. 為了保證IVIM參數的一致性采用與NLLS方法一致的擬合邊界. 式(5)中干擾參數S(b0)和將通過一個先驗分布p(S(b0)來消除. 先驗分布如式(6)定義所示:

2 材料和方法

2.1 數據集及數據處理

本文真實大腦測試臨床數據集采集于哈爾濱腫瘤醫(yī)院, 共計5名受試者. 所有磁共振圖像均采集于飛利浦3TMRI掃描儀(Philips MR 53.0, 3.0 Tesla MRI scanner), 采集序列為平面成像自旋回波序列(EPI SE).其采集參數如下: 分辨率(matrix size)為128×128、重復采集時間(TR)為3 000 ms、回波時間(TE)為70.28 ms、切片厚度(slice thickness)為2.5 mm. 擴散敏感度b值序列為(0, 3, 5, 10, 25, 50, 75, 100, 200, 400, 600,800 s/mm2). 本文首先利用MRIcron軟件對原采集數據進行格式轉換, 再使用Local PCA算法對數據切片上進行降噪. 隨后使用SyN算法將各個數據配準到b0圖像上[31], 最后為了對大腦不同感興趣區(qū)域進行分析, 使用FSL工具對受試者大腦進行去頭骨處理[32], 再使用SPM8對灰質、白質、腦脊液進行了分割.

2.2 蒙特卡洛模擬仿真

為了訓練和對比傳統機器學習方法、ANN以及DCNN方法對于IVIM參數預測的性能, 我們基于蒙特卡羅模擬(Monte-Carlo simulation)方法, 分別仿真了訓練模擬信號和測試仿真數據. 訓練模擬信號是在真實大腦IVIM參數臨床表現的限定范圍里, 均勻隨機選取各IVIM參數的標簽再利用式(1)和與真實數據集相同的擴散敏感度b值序列, 得到模擬信號數據其中模擬擴散參數取值范圍為(0.5, 3×10-3mm2/s), 模擬灌注參數取值范圍為(0.003, 0.1 mm2/s), 模擬灌注體積分數為(0, 1). 隨后, 再向每個模擬信號中加入萊斯噪聲(Rician noise)模擬出與真實數據相同分布情況的訓練仿真信號添加噪聲方法如式(7)所示:

其中,GN1和GN2滿足均值為0、方差為1 /SNR的兩不相關高斯噪聲. 為了防止網絡模型過擬合, 增加訓練模擬信號的復雜性和多樣性,SNR在(10, 120)區(qū)間隨機取值.

由于IVIM數據缺乏參數的金標準, 參數估計的準確性不能直接從真實數據中確定, 但可以通過計算機模擬出測試仿真數據來評估IVIM參數估計的準確性.為了得到真實的IVIM參數標簽來評估各方法的性能,我們將一例真實數據的大腦, 利用SPM8軟件分割得到了正常人腦組織的1個灰質(gray matter, GM)、2個白質(white matter, GM)和1個腦脊液(cerebrospinal fluid, CSF)共計4個區(qū)域.

將各個區(qū)域賦予設定的真實參數, 如下所示(灰質(0.002, 0.07, 0.4), 白質區(qū)域1 (0.001, 0.02, 0.06), 白質區(qū)域2 (0.001, 0.03,0.08), 腦脊液(0.002, 0.07, 0.4), 每種特定組織的參數值可以涵蓋文獻和研究中報道的典型值范圍. 測試仿真數據的IVIM參數圖分別如圖1各子圖所示, 同時利用式(7)向信號中添加噪聲, 信噪比為(20, 30, 40, 50,60).

圖1 測試仿真參數圖

2.3 網絡結構

近年來, 利用人工神經網絡(artificial neural network,ANN)對IVIM數據進行分析已經被報道. 然而, 這種方法并沒有完全考慮IVIM參數與b值之間的關系, 也沒有體素信號間上下文信息對IVIM參數的影響. 針對于此問題, 根據IVIM磁共振信號的特點, 利用動態(tài)卷積模塊(dynamic convolution module, DCM)和殘差網絡(residual network)相結合的方式, 以提高IVIM參數的估計準確度和良好的視覺效果.

在傳統的卷積中, 每個神經元的輸入被描述為y=σ(ωx+b), 其中,ω代表著卷積核,σ是激活函數,b是神經元的偏置. 動態(tài)卷積核是由一系列的并行卷積核ωi和 相應的線性權重ηi聚合而成, 它能夠在不增加網絡深度的情況下, 通過學習于通道的線性權重ηi, 也就是根據IVIM信號衰減特性, 自適應調整卷積核參數, 增加卷積神經網絡模型的表達能力. 動態(tài)卷積核聚合過程如式(8)所示:其中,N表示線性權重 ηi和 并行卷積核ωi的個數. 動態(tài)卷積包含權重預測模塊(weight prediction module, WPM)和卷積聚合模塊(convolution aggregation module,CAM). WPM系數預測模塊由全局平均池化層(global average pooling, GAP)、全連接層、激活函數和Softmax組成. IVIM數據首先通過全局平均池化, 獲得了由于b值導致信號衰減的整體信息, 再通過全連接層和Softmax得到線性權重 ηi. Softmax的運用能夠保證線性權重ηi滿足式(8)的條件.

此時得到的線性權重ηi, 通過卷積聚合模塊與對應的并行卷積核ωi進行乘積再相加來實現, 如式(8)所示聚合過程, 動態(tài)卷積流程圖如圖2所示.

圖2 動態(tài)卷積模塊流程圖

上述動態(tài)卷積的性質表明, 將原本固定且靜態(tài)的卷積核替換為動態(tài)卷積核, 能夠增強網絡對于b值關于IVIM參數貢獻程度的表達能力, 同時使用一維卷積神經網絡能夠學習到相鄰信號中的上下文信息. 因此,我們設計了一個基于動態(tài)卷積模塊的一維卷積神經網絡(dynamic convolutional neural network, DCNN)來估計IVIM參數, 詳細的結構如圖3所示. DCNN由6組一維動態(tài)卷積模塊組成, 考慮到網絡層數的加深會導致梯度消失的問題, 其中3組分別使用了殘差網絡. 殘差網絡可以將輸入通道的IVIM信號衰減信息傳遞到輸出通道, 使得殘差網絡中間的卷積層只需要學習殘差, 避免信息遺漏.

圖3 基于動態(tài)卷積的IVIM參數估計網絡圖

在DCNN網絡中, 利用初始信號S(b0)對IVIM擴散加權信號S(bn)進行歸一化操作得到S(bn)/S(b0), 以格式為(B,N,L)的數據作為網絡的輸入, 其中,B是batch_size,N是b值個數,L是體素個數. 在動態(tài)卷積模塊中, 網絡分別采用了核大小為3的一維動態(tài)卷積和1的一維常規(guī)卷積作為IVIM參數擬合主體. 將歸一化后的信號S(bn)/S(b0)輸入到核大小為3的一維動態(tài)卷積模塊中, 能夠捕獲相鄰體素信號間上下文信息, 然后輸入到核大小為1的一維常規(guī)卷積中, 以提高線性組合的擬合能力. 然后, 數據通過一個帶有Sigmoid激活函數的動態(tài)卷積的網絡輸出層, 其卷積核為3, 輸出通道數為3. 最后, 將這3個輸出向量重塑為3張圖像,分別代表預測的3個IVIM參數. 另外, 輸出神經元的參數用式(9)進行歸一化, 以確保估計的IVIM參數在合理范圍內.

動態(tài)卷積神經網絡(DCNN)將通過最小化預測信號和真實的歸一化信號x之間的均方差來學習網絡參數, 如式(10)所示:

3 實驗設計

3.1 實驗環(huán)境及訓練細節(jié)

本實驗采用的硬件環(huán)境, 數據預處理部分在顯卡為NVIDIA-GTX 1050ti, 顯存共4 GB, CPU為Intel(R)Core(TM) i7-8700, 內存共36 GB的主機上, 使用PyCharm作為本實驗項目管理工具. 深度學習網絡在擁有32 GB的NVIDIA Tesla V100 GPU服務器上, 使用Python 3.6語言、PyTorch 1.5工具包及Cuda 10.1驅動共同支撐的軟件環(huán)境下進行加速訓練. 本實驗所有網絡模型, 均使用學習率為1×10-3的Adam作為優(yōu)化器, 其中batch_size為40, 共訓練2 000個Epoch. 為了防止過擬合, 網絡訓練使用了早停機制(early stopping),當連續(xù)15個訓練周期未下降則停止訓練. 為了比較和量化提出方法的可重復性和穩(wěn)定性, 使用相同的訓練仿真數據進行了10次訓練.

3.2 評價指標

為了評價不同參數估計方法對IVIM參數估計的性能, 將采用平均值(Mean)、標準差(standard deviation,Std)、相對誤差(relative bias,RB)、變異系數(coefficient of variation,COV)和均方根誤差(root mean square error,RMSE)作為定量標準. 其中, 由于真實圖像沒有絕對準確的真值標簽, 我們將用平均值±標準差(Mean±Std)和COV變異系數來評定真實圖像IVIM參數. 對于仿真圖像而言, 我們計算出仿真圖像的RB相對誤差、COV變異系數和RMSE均方根誤差.

其中, σt代表著不同信噪比下各參數圖的方差,S是添加噪聲圖像的區(qū)域個數,T代表著感興趣區(qū)域的體素點個數. 其中,COV變異系數能夠反映在相同結構下的參數均勻程度, 理想的COV為0. 高COV值表示有高噪聲的映射, 圖像視覺質量低.

4 實驗結果及分析

為了評估本文模型的有效性和穩(wěn)定性, 本文與NLLS、Bayes、ANN三種不同方法對仿真測試數據集和真實的采集數據集進行了IVIM參數估計, 實驗可視化和定量結果如下.

4.1 仿真測試數據

圖4顯示了4種IVIM估計方法對測試仿真數據在SNR為20, 30, 40, 50, 60上的參數圖結果. 很明顯, 基于動態(tài)卷積的方法大大提高了3個IVIM參數圖的圖像質量. 此外, 與傳統機器學習方法(NLLS、Bayes)相比, DCNN更有效地抑制了噪聲的生成, 尤其在灌注相關參數圖(Dp,fp)上, 有更清晰的紋理信息, 空間信息也損失得更少, 同時在均勻區(qū)域更加平滑. 圖5是各方法對于IVIM測試仿真參數圖在不同SNR下的評價指標結果圖. 由圖5可知, DCNN在Dp參數圖的估計上全面優(yōu)于NLLS、Bayes和ANN, 說明該方法在視覺質量的優(yōu)越性和數值上的合理性. 對于Dp和fp而言,RMSE和COV有大幅度的下降, 說明采用一維動態(tài)卷積神經網絡在考慮上下文信息有較強的特征提取能力, 使得體素與體素點之間有更強的相關性和更強大的抗噪聲能力. DCNN的RB大幅下降, 說明動態(tài)卷積模塊能夠學習到來自于b值的貢獻權重信息, 使得不再像NLLS那樣高估灌注系數數值. 除了NLLS之外, 各方法對于Dt參數圖的評價指標相差不大. 但值得注意的是, 在SNR大于40的時候, Bayes的RB比DCNN小約0.25%.

圖4 各方法的測試仿真數據參數結果圖

圖5 各方法IVIM測試仿真參數結果在不同SNR下的評價指標圖

4.2 真實采集數據

圖6為一健康受試者大腦的4種IVIM參數估計算法的可視化結果. 與預期的一樣, 用不同方法估計的IVIM參數的變化與仿真結果一致. 整體來說, DCNN提供了最平滑、噪聲最小的參數映射. 具體而言, 各種方法對于Dt的參數估計最為相似, 但是NLLS的Dt參數圖顯然受到了強烈的噪聲影響. 兩種深度學習方法(ANN, DCNN)的Dt的參數圖相對于傳統機器學習方法(NLLS, Bayes)相對干凈、平滑, 但DCNN在灰質上的紋理結構有著更好的視覺效果. 對于灌注相關系數(Dp和fp)而言, 不同估計方法差異較為明顯. 對于Dp而言, 可以看出NLLS存在著大量的異常值. 盡管貝葉斯方法改善了Dt和fp參數映射, 但在Dp映射中仍有許多異常值. 這是因為Bayes將NLLS估計結果作為先驗, 而這種由于先驗的影響導致仍然存在許多離群值是不可避免的. 對于DCNN而言, 所以生成的Dp圖受異常值的影響較小, 視覺上的紋理信息豐富,明顯優(yōu)于ANN神經網絡. 綜上所述, 我們發(fā)現DCNN可以更有效地減少異常值, 并且可以保留大部分的空間和紋理信息, 以獲得最佳的視覺質量. 這驗證了所提出的基于DCNN的方法可以提高IVIM參數映射的視覺質量. 除可視化結果之外, 表1總結了用不同擬合方法估計的WM和GM的IVIM參數的平均值±標準差和COV變異系數. 對于WM和GM, 不同擬合方法的參數圖Dt的平均值和Std之間沒有顯著差異. 然而, 由所提出的動態(tài)卷積方法得到的參數映射Dp和fp的Std和COV比其他方法要小得多, 說明該方法具有IVIM參數估計的最優(yōu)穩(wěn)定性.

表1 各方法的IVIM參數定量分析結果

圖6 各方法對IVIM參數結果的可視化圖

4.3 一致性與收斂性分析

本文利用相同的訓練仿真數據集, 采用隨機初始化對兩種深度學習方法(ANN和DCNN)獨立訓練了10次, 以研究不同深度學習方法性能的一致性. 圖7箱型圖顯示了5名健康受試者在腦灰質和白質上的COV變異系數分布情況. 總的來說, DCNN和ANN的方法產生的Dt參數圖有著相似的數據分布情況, 無論是在白質還是在灰質. 雖然通過DCNN模型得到的COV均值偏移量略高于ANN, 但在灰質上的四分位數范圍小于ANN. 然而, 相比于ANN, DCNN得到的灌注相關參數(Dp和fp)的COV分布有顯著差異. 與仿真測試集和真實采集數據結果類似,Dp參數圖在視覺上展現不同的分布. 在Dp上, ANN獲得的平均COV和四分位數范圍明顯高于DCNN, WM和GM獲得的平均COV均超過75%, 而DCNN的平均COV, GM為45.9%,WM為26.3%. 對于fp而言, ANN獲得的平均COV在GM上為88.3%, WM上為70.8%, 而DCNN在GM上為76.9%, WM上為53.1%. 實驗結果表明, 基于動態(tài)卷積DCNN的方法的性能比ANN具有更好的一致性,特別是對灌注相關參數的估計, 這表明我們提出的方法更加穩(wěn)定.

圖7 深度學習算法一致性箱線圖

圖8顯示了不同深度學習網絡模型前500次訓練損失的收斂曲線圖. 由圖8可知, DCNN大約在130次迭代后, 目標優(yōu)化函數數值趨于穩(wěn)定, 約為0.004 0. 然而, ANN收斂曲線有較大幅度的振蕩現象, 在300次迭代后損失函數數值約為0.013 4, 說明ANN方法陷入局部解且逃逸能力不佳. 上述表明DCNN模型的收斂性和優(yōu)化目標損失函數的效果要好于ANN方法.

圖8 網絡模型訓練損失收斂曲線圖

5 結論

在本研究中, 我們提出了一個基于動態(tài)卷積的模型來提高IVIM參數估計的精度, 該模型同時考慮了不同b值的DW信號對IVIM參數估計的貢獻程度以及相鄰體素之間的上下文信息. 在模擬仿真和真實采集的數據集上, 實驗結果表明, 相比于NLLS、Bayes和ANN方法, 本文模型估計IVIM參數的能力更佳. 該方法可以生成具有更小COV變異系數和更少離群值的IVIM參數圖, 尤其是對于Dp參數圖. 此外, 本文還考慮到方法的一致性, 重復實驗結果表明, 本文所提出的方法比傳統ANN神經網絡估計結果更穩(wěn)定, 滿足了臨床應用的實際需要. 接下來的工作中, 我們將會在更多不同疾病的受試者來驗證我們的網絡模型的性能評估.同時, 我們將會考慮訓練一個泛化網絡, 可以接受任何器官的任何采集協議獲得的DW信號, 進而實現不同采集協議下IVIM 參數的準確估計.