基于深度先驗及非局部相似性的壓縮感知核磁共振成像

宗春梅，張月琴，曹建芳，趙青杉

（1.忻州師范學院計算機系，山西忻州 034000；2.太原理工大學計算機科學與技術學院，太原 030024）

（*通信作者電子郵箱zongcm123@163.com）

0 引言

壓縮感知核磁共振成像（Compressed Sensing Magnetic Resonance Imaging，CSMRI）是指通過觀測的高度欠采樣的k空間數據（圖像的傅里葉變換，也稱測量數據）重建原始圖像［1］。在采樣矩陣滿足某種條件下，該技術能夠從遠低于奈奎斯特采樣率的k空間數據中精確恢復原始圖像［2］。CSMRI理論一經提出就以采樣時間短、存儲數據少的優點得到了學者們的關注。如何從少量的測量數據中高質量地恢復原始圖像是CSMRI 領域面臨的一大挑戰。CSMRI 僅利用部分測量數據即包含少量信息的觀測數據進行圖像重建，該問題有無數個解，是一個病態問題。解決該問題的一個有效方式是利用圖像固有的先驗知識進行圖像重建。近些年來，國內國外涌現出了大量的CSMRI 算法。根據算法利用的先驗知識大致可以分為三大類：基于稀疏性的算法［3-5］、基于非局部相似性的算法［6-8］、基于深度先驗的算法［9-12］。

基于稀疏性的算法是利用圖像在變換域或者字典下的稀疏性進行圖像重建。其中，稀疏性又可以分為全局稀疏性和局部稀疏性。全局稀疏性是指圖像在梯度域、小波變換域、緊標架下的稀疏性，而局部稀疏性是指圖像塊在固定字典或自適應字典下的稀疏性。基于全局稀疏性的方法重建速度相對較快，但是重建質量有待提高。基于局部稀疏性的方法通過字典學習技術進行圖像重建［5］，利用交替優化的方式在迭代過程中聯合優化圖像和字典。相較于基于全局稀疏性的方法，該類方法能夠學習自適應字典，重建質量高，但成像速度相對較低。

基于非局部相似性的CSMRI 算法是利用非局部區域內圖像塊之間的相似性進行重建，本質是圖像的當前塊估計通過一定區域內的相似塊加權得到。近些年，利用塊匹配三維濾波（Block Matching and 3D Filtering，BM3D）高斯去噪器［13］隱式利用非局部相似性進行圖像重建得到了學者們的關注。Eksioglu［6］利用BM3D 去噪器構建了去耦合的CSMRI 框架，該框架能夠隱式地利用BM3D 去噪器包含的圖像先驗知識，即圖像的非局部相似性及圖像在三維小波變換下的稀疏性。Shi 等［7］利用BM3D 去噪器構建了正則化模型，并用于提升基于局部稀疏性的CSMRI 算法重建質量。實驗結果表明該類方法能夠獲得較高質量的圖像。

基于深度先驗的算法是通過深度神經網絡進行圖像重建。Yang 等［9］利用交替方向乘子法（Alternating Direction Method of Multipliers，ADMM）求解構建的CSMRI 優化問題，并將求解框架展開成網絡，提出了ADMM-CSNet（Alternating Direction Method of Multipliers-Compressed Sensing Network）算法。ADMM-CSNet算法能夠通過事先訓練好的深度神經網絡進行快速成像且取得較好的重建效果。基于深度先驗的算法雖然能夠取得較好的性能，但是該類直接訓練端到端深度神經網絡的方法泛化能力差，且不具可解釋性。

上述算法雖然取得了較好的重建效果，但在低采樣率時，重構質量仍然有提升空間。不同于上述算法，本文利用深度先驗和非局部相似性進行壓縮感知核磁共振成像。本文利用深度去噪器和BM3D 去噪器通過即插即用模型將這兩種先驗知識引入到圖像重建中。受基于深度先驗的稀疏表示模型［14］啟發，本文提出基于兩種去噪器的稀疏表示模型。該模型能夠融合深度先驗、圖像的非局部相似性及圖像在三維變換下的稀疏性進行圖像重建。結合該模型與數據保真模型構建CSMRI 優化問題，利用交替優化的方式對該問題進行有效求解。本文利用多種先驗知識進行成像期望得到高質量的圖像重建。

1 基于稀疏性的壓縮感知核磁共振成像

假設原始核磁共振圖像為x∈RN，CSMRI 采樣模型可描述為

其中：Fu∈Cm×N（m＜N）表示欠采樣傅里葉編碼矩陣，也即欠采樣算子；ratio=m/N定義為采樣率；y∈Cm表示觀測的k空間數據即測量數據；n∈Cm表示高斯噪聲向量。如何利用測量數據y重建原始圖像x是CSMRI算法所要解決的問題。

基于全局稀疏性的CSMRI算法通常構建以下優化問題：

其中：ψ是某種全局變換；λ為正則化參數。||·||p表示lp范數：p=0時，即為l0偽范數，衡量對應向量中非零元素的個數；當p=1 時，||·||p為l1范數，表示對應向量中每個元素絕對值的和。上述優化問題可通過梯度下降法、變量分離法或其加速方法進行求解。基于全局稀疏性的方法能夠快速進行圖像重建，但由于全局變換是固定的，重建質量有待提升。

不同于上述算法，基于局部稀疏性的CSMRI 算法通常利用自適應字典進行重建。在重建過程中，通過估計圖像學習自適應字典。該類方法通常構建以下優化問題［5］：

其中：?是圖像塊索引的集合；Ri∈RB×N表示取塊算子；γ是平衡參數；D∈CB×K表示待訓練的塊字典；αi∈CK表示圖像塊Rix在字典D下的表示系數；Γ是所有表示系數的集合；T為稀疏度。上述優化問題通常利用交替優化的方式分字典學習和圖像更新兩步進行求解。由于字典能夠通過測量數據進行訓練，學習的字典能夠對重建圖像自適應稀疏表示。基于字典學習的CSMRI 算法較基于全局變換的CSMRI 算法重建質量高。

為進一步提升重建質量，Eksioglu［6］利用BM3D 去噪器隱式利用隱含在BM3D 去噪器中的先驗知識進行圖像重建，提出了基于三維塊匹配的核磁共振成像算法BM3D-MRI（Block Matching and 3D filtering-Magnetic Resonance Imaging）算法。該算法構建以下去耦合的CSMRI 優化框架（對于第t次迭代）：

其中：ΦBM3D表示BM3D 去噪器對應的解析稀疏算子［6］。式（4）可通過求導并令導數為0 獲得閉式解，式（5）直接利用BM3D 去噪器近似求解。由于上述框架在優化過程中利用BM3D 去噪器對估計圖像進行去噪處理，該框架能夠隱含地利用非局部相似性。

受上述算法及基于深度先驗的稀疏表示模型［14］啟發，本文通過BM3D 去噪器和深度去噪器相結合，提出融合多種先驗知識的CSMRI 算法。本文構建了基于多種去噪器的稀疏表示模型，該模型能夠將多種去噪器蘊含的先驗知識融合到圖像重建中。多種先驗知識的利用能夠提升低采樣率下的重建質量。

2 基于多種去噪器的稀疏表示模型

如何利用多種先驗知識進行CSMRI 是本文解決的一個關鍵問題。本文利用兩種不同的去噪器構建稀疏表示模型以利用多種先驗知識。傳統稀疏表示模型假定圖像x在標架或某種變換W下是稀疏的，通常構建以下優化模型以求解稀疏表示系數［15］：

其中參數η控制稀疏程度。對于任意變換W，上述優化模型可以通過迭代閾值方法求解。為了降低計算復雜度，本文將緊約束施加到標架W上，即W為緊標架WTW=I。基于緊標架的優化問題（6）可直接通過簡單的閾值方法求解。圖像經過BM3D 去噪器的濾波圖像在緊標架下的表示系數與原始圖像在同一緊標架下的表示系數應該是近似的，兩個系數向量的差值應是稀疏的。基于該認識，本文構建基于多種去噪器的稀疏表示模型，優化模型可表示為：

其中：De(x)表示利用去噪器對圖像x進行濾波；μ表示稀疏正則化參數，用以控制濾波圖像在緊標架下的稀疏程度。De(x)可以是單個去噪器，也可以是多種去噪器的加權求和即表示第j個去噪器對圖像進行濾波以后的圖像）。文獻［14］指出利用不同的高斯去噪器可以將不同的先驗知識引入到圖像重建中。本文利用兩種去噪器構建該模型，即深度去噪器與BM3D 去噪器。在這種情況下，濾波圖像可表示為。其中，Deep(x)表示利用深度去噪器對圖像進行濾波，BM3D(x)表示利用BM3D 去噪器對圖像進行濾波。本文聯合采用這兩種去噪器，試圖融合兩種去噪器在去噪或濾波過程中利用的先驗知識進行圖像重建。深度高斯去噪器利用深度先驗進行圖像去噪，利用該去噪以將深度先驗引入到圖像重建中。深度去噪器通常未充分利用圖像的自相似性，BM3D 去噪器彌補了該不足，該去噪器通過利用非局部區域內的相似性先驗知識進行圖像去噪。本文利用這兩種去噪器進行聯合去噪或濾波以期望利用深度先驗和非局部相似性先驗的互補性進行高質量圖像重建。

式（8）可通過對濾波圖像在緊標架下的表示系數進行閾值處理求解。如果p=1，則，此時表示閾值為的軟閾值算子，定義為T[·，ε]=soft{·，ε}=sign(·)max(|·|-ε，0)，sign(·)為符號函數，max(·)為取最大值算子。如果p=0，則，此時表示硬閾值算子，定義為：如果，否則。需要指出的是，閾值算子操作是逐元素進行的，因此上述對向量的處理是對每一元素進行閾值處理。

構建的上述正則化模型R(x)是一個包絡函數，可以直接作為正則項引入到成像的優化模型中。式（9）中的第一項表示圖像在緊標架下的表示系數近似為α，第二項表示該系數逼近于輔助系數向量。其中該輔助系數向量的更新能夠通過兩種去噪器將深度先驗、非局部相似性及三維變換下的稀疏性隱式地進行利用。由于式（9）中的輔助系數向量可通過求解式（8）獲得，故基于BM3D 去噪器的稀疏表示正則化模型可寫成以下簡潔的形式：

構建的正則化模型（10）具有以下優點：

1）有效。該正則化模型能夠通過多種去噪器利用多種互補先驗知識進行圖像重建，多種先驗知識有益于圖像重建；

2）靈活。該模型的靈活性主要體現在以下兩個方面：第一，去噪器可以換成其他有效的去噪器以利用不同的先驗知識進行圖像重建；第二，該正則化模型可用于其他成像應用，不局限于壓縮感知核磁共振成像。

3 融合深度先驗及非局部相似性的CSMRI

本文利用基于兩種去噪器的稀疏表示模型進行壓縮感知核磁共振成像以解決現有CSMRI 算法重建質量低的問題，構建以下優化模型：

式（11）代價函數中：第一項為數據保真項，保證重建的圖像與測量數據y相匹配；第二、三項為本文構建的基于多種去噪器的稀疏表示模型。由于上述代價函數采用l2范數和l1范數，如果固定輔助系數向量T[WDe(x)，ε]，所對應的優化是一個凸優化問題。然而，由于閾值算子（閾值算子是非凸的）和多種去噪器（去噪器可看成去噪函數，通常是非凸的）存在于最后一項，上述優化問題在考慮x和α聯合優化時，該問題是一個復雜的非凸優化問題。為了對其進行有效求解，本文假設濾波圖像De(x)近似于估計圖像的濾波圖像De(x(t-1))，即De(x)≈De(x(t-1))。基于該假設，優化問題（11）可改寫為（對于第t次迭代）：

利用交替優化的方法求解上述優化問題，對于第t次迭代，分為以下兩步對問題（12）進行求解：

1）系數更新步驟，固定估計圖像x(t-1)，更新稀疏系數α(t)的子問題為：

上述優化問題可通過軟閾值算子進行求解：

實驗結果表明，T[·，ε1]采用硬閾值算子成像效果優于利用軟閾值算子，因此本文采用硬閾值算子。對濾波圖像在緊標架下的表示系數進行硬閾值處理，得到輔助系數向量T[WDe(x(t-1))，ε1]。

2）圖像更新步驟，固定稀疏系數α(t)，更新圖像x(t)的子問題為：

對式（15）中的代價函數求導，并令導數為零：

傅里葉變換矩陣F∈CN×N通常為酉矩陣，滿足FHF=I。根據這一等式，式（16）可改寫為：

其中：Fx(t)(kx，ky)表示在位置(kx，ky)的更新值；Ω是被采樣的k空間數據的集合。最后對Fx(t)進行傅里葉反變換，獲得最終圖像重建結果。

交替優化稀疏系數與圖像直到達到終止條件可以得到優化問題（11）的一個近似解。本文算法如下：

4 實驗與結果分析

為了驗證本文算法的有效性，本文采用了6 幅大小為256 × 256 的核磁共振圖像作為測試圖像進行壓縮感知核磁共振成像實驗，測試圖像分別為Brain、Shoulder、Bone、Head、Bust、Brain2。采用偽隨機采樣算子對k空間數據進行欠采樣，圖1給出了采樣算子和6幅原始圖像。所有算法均在配置為Core i7-7700 主頻3.6 GHz CPU，內存8 GB 的PC 上進行測試，軟件平臺為Windows 10 64 位操作系統，Matlab 2018b。為了模擬更真實的采樣環境，本文對采樣值施加噪聲標準差為σ=1 的噪聲。本文算法在不同采樣率情況下與基于小波樹稀疏性的核磁共振成像算法WaTMRI（Magnetic Resonance Imaging with Wavelet Tree sparsity）［4］、基于字典學習的核磁共振成像算法 DLMRI（Dictionary Learning for Magnetic Resonance Imaging）［5］、基于字典更新及塊匹配和三維濾波的核磁共振成像算法（Magnetic Resonance Imaging based on Dictionary Updating and Block Matching and 3D filtering，DUMRI-BM3D）［7］、BM3D-MRI 算法［6］、廣義近似消息傳遞去噪算法DAMP（Denoising Approximate Message Passing）［8］進行對比，對比算法均采用軟件包中給出的默認參數。本文分別從客觀角度和主觀視覺角度對測試算法進行對比。采用峰值信噪比（Peak Signal-to-Noise Ratio，PSNR）作為評價重建質量好壞的客觀標準，主觀視覺評價則計算特征相似性（Feature SIMilarity，FSIM）［16］，同時給出各個算法的部分重建圖像進行對比。PSNR 值和FSIM 值越大表明重建圖像視覺效果越好，重建質量越高。

在實驗中，本文的深度去噪器選用近些年比較優秀的深度去噪器FFDNet（Fast and Flexible Denoising convolutional neural Network）［17］。輸入去噪器的噪聲標準差采用BM3DMRI 算法給出的非線性構建方法［6］。假設該方法計算的標準差為σ，為實現較好的重建效果，輸入BM3D與FFDNet去噪器的噪聲標準差分別為6σ、3σ。正則化參數λ、閾值ε1及ε2均采用PSNR 值最大的準則進行調節，對于某一幅測試圖像，固定一個參數，去調節另一個參數使得獲得的PSNR 值最大。當調到最大效果時，固定該參數，去調節另一個。閾值ε1與ε2應與估計圖像中的噪聲標準差有關，經驗表明，當λ=0.01、ε1=0.5σ、ε2=5σ時，圖像重建效果最好。當參數調整好以后，利用這套參數對不同采樣率下、不同圖像進行核磁共振成像。緊標架選用離散余弦變換并使用文獻［14］公布的默認參數。為了公平比較，提出的算法對于不同的采樣率選用同一套參數。為了說明多種互補先驗的有效性，本文將僅利用FFDNet 去噪器的方法作為基準算法，簡記為CSMRIFFDNet。

圖1 采樣算子和原始圖像Fig.1 Sampling operator and original images

表1 給出了在采樣率為0.02、0.06、0.09 及0.13 的情況下測試算法的PSNR 值（由于版面限制，本文僅給出5 種具有代表性的測試算法結果），從中可以看出：1）對于多數情況，本文算法的PSNR 值最高；2）以Brain 圖像為例，在采樣率為0.02 情況下，本文算法重構圖像的PSNR 值比DUMRI-BM3D算法、BM3D-MRI 算法、DAMP 算法及CSMRI-FFDNet 算法分別提高了1.05 dB、0.45 dB、0.85 dB 及0.34 dB；3）在低采樣率下，本文算法仍然能夠重構出較高PSNR 的圖像；4）利用兩種去噪器的本文算法優于僅利用FFDNet 去噪器的CSMRIFFDNet 算法。CSMRI-FFDNet 算法利用了深度先驗，重建平均PSNR 值優于其他對比算法。深度先驗并未考慮圖像非局部區域的相似性，由于本文算法能夠利用深度先驗與非局部相似性進行圖像重建，因此本文算法優于僅利用深度先驗的CSMRI-FFDNet算法。

為了進一步衡量本文算法主觀視覺的優點，表1 給出在不同采樣率下測試算法獲得的FSIM 值。從表1 中可以看出，本文算法對于不同圖像、不同采樣率都能夠獲得最高的FSIM值，也就意味著視覺效果最好。圖2 給出了測試算法對圖像Head的重建結果，從重建圖像中可以看出：WaTMRI的重建質量最差、最模糊，含有大量的噪聲；DLMRI 算法重建圖像消除了部分噪聲，但仍然丟失了大量細節信息，并且重建圖像具有塊效應；DUMRI-BM3D 算法、BM3D-MRI 算法及DAMP 算法消除了噪聲，但是重建圖像丟失了部分細節信息，具有明顯的偽影效應；而本文算法消除了偽影，保留了大量細節信息，視覺效果最好。

表1 不同算法PSNR值和FSIM值比較Tab.1 Comparison of PSNR and FSIM values of different algorithms

為了衡量算法的成像速度，表2 給出了測試算法的平均運行時間，是測試算法對4 種采樣率、6 幅圖像的運行時間平均值。從表2 可以看出，僅利用全局變換的WaTMRI 算法耗時最少，利用字典學習進行成像的DLMRI 算法用時最多。DUMRI-BM3D 算法利用字典更新和BM3D 去噪器進行成像，因此成像時間相對較長。本文算法的運行時間較WaTMRI算法、BM3D-MRI算法及CSMRI-FFDNet長。這是因為本文算法利用了稀疏表示和BM3D 去噪器。為提升成像速度，本文算法可通過并行計算進行加速。雖然本文算法比這3 種算法耗時長，但本文算法重建質量高于這3 種算法，更加適用于低采樣率下對重建質量要求高的場合。

圖2 測試算法的Head重建圖像（采樣率為0.06）Fig.2 Reconstructed head images by different test algorithms（sampling ratio of 0.06）

表2 不同算法的運行時間比較Tab.2 Running time comparison of different algorithms

5 結語

本文針對深度先驗和非局部相似性的互補性，利用基于多種去噪器的稀疏表示模型，提出了融合兩種先驗知識的壓縮感知核磁共振圖像重構算法。該算法將基于多種去噪器的稀疏表示模型結合到CSMRI 的代價函數中，利用交替優化方法有效地求解了所對應的優化問題。實驗結果表明，提出的算法能夠在低采樣率下獲得較高的重建質量，且兩種先驗知識的融合能夠提升圖像的重建質量。如何加速本文算法是以后的研究方向。