基于xy非線性梯度場磁共振成像研究

徐雅潔， 楊曉冬

(1.中國科學院 蘇州生物醫(yī)學工程技術研究所， 江蘇，蘇州 215163；2. 中國科學院大學， 北京 100049)

?

基于xy非線性梯度場磁共振成像研究

徐雅潔1，2， 楊曉冬1

(1.中國科學院 蘇州生物醫(yī)學工程技術研究所， 江蘇，蘇州 215163；2. 中國科學院大學， 北京 100049)

為解決非線性梯度成像中心分辨率低的問題，提出了一種基于xy非線性梯度場的成像方法，并討論了非線性梯度場中心分布曲線對磁共振圖像的影響. 研究結果表明，在半徑r=1.5時，圓形中心分布曲線的圖像重建結果最優(yōu)；與O-Space技術對比實驗表明，xy非線性梯度場在不同加速因子R下都得到了更好的重建結果.

非線性梯度場成像；中心分布曲線；CG迭代算法

非線性梯度場成像是弗賴堡大學的Juergen Hennig教授在2008年提出的磁共振成像新技術，以局部梯度場并行成像技術[1](parallel imaging technique using local gradients，Patloc)為代表，利用非線性空間編碼磁場(spatial encoding magnetic field, SEM)和多通道并行采集方法實現磁共振并行成像. Patloc中的曲面SEMs更適應被試結構，在提高編碼效率、縮短成像時間、提供更多的診斷學信息等方面具有很大潛力. 隨著高場磁共振設備的普及，梯度線圈中的感應電流對人體神經末梢產生的電刺激導致的被試不舒適感也逐漸增加，稱為外周神經刺激(peripheral nerve stimulation，PNS). 相關研究表明，人體所承受PNS與磁場強度B和磁場變化率dB/dt都有一定的聯系[2]，其中前者更為密切，而非線性梯度場中磁場強度的絕對值比線性梯度場小，因而在降低被試PNS方面也具有較大優(yōu)勢.

然而Patloc技術是使用雙曲線型磁場x2-y2和反比例函數型磁場xy(兩者僅存在相位差)，采取傳統的K空間編碼方式進行的成像. 由于非線性梯度場具有邊緣編碼梯度高，中心梯度低的特性，而該方法采用的陣列接收線圈為成像區(qū)域外周對稱排列式，使得其在成像中心的靈敏度也相對較低，造成采用該方法得到的成像結果通常邊緣分辨率相對較高成像清晰，而中心卻比較模糊. 為了解決這個問題，研究者們嘗試利用多SEMs進行空間編碼. Gallichan[3]設計的4DRIO(4-dimensional radial in/out)技術使用兩組空間編碼磁場x，y和x2-y2，xy，采用輻射狀采樣進行成像，明顯改善了中心部分分辨率，但輻射角之間的關系對成像結果影響較大；隨后，Stockmann[4-6]提出了O-Space 技術，利用Z2階磁場分量：2z2-(x2+y2)，通過線性梯度場改變Z2磁場中心的位置，在整個成像空間均勻覆蓋編碼磁場，從而解決了非線性梯度場中心模糊的缺陷.

從Z2函數表達式可知，此非線性磁場隨三坐標軸x，y，z變化，在3D成像中由于z2場變化會導致散相[4]，因而對成像存在較大影響.x2-y2和xy非線性梯度場僅與x，y坐標軸有關，故應用此磁場比O-space技術中的Z2 階磁場更為方便，且O-Space技術根據其磁場結構可知僅能應用在橫斷面成像中，對冠狀面，矢狀面沒有相對應的球諧函數磁場分量，而反比例函數型梯度場則在x，y，z方向都有可利用的非線性梯度場，分別為zy，zx，xy. 本文即以xy磁場為代表研究非線性梯度場成像技術.

1 電磁場理論

根據Maxwell方程組[7]，通過對標量磁位的變換求解和球諧函數展開，有

(1)

笛卡爾坐標系下展開式為

(2)

圖1中給出了2z2-(x2+y2)，zy，zx，xy，x2-y2磁場示意圖[8]，zy，zx，xy磁場分別關于x，y，z方向對稱，xy與x2-y2階磁場僅存在45°相位差. 在本文應用中zy，zx，xy磁場可以分別實現冠狀面、矢狀面、橫斷面的成像，下面以xy磁場為例介紹非線性梯度場成像.

陣列線圈接收信號可以寫成[9]

(3)

式中：α為并行接收線圈；m(x)為感興趣區(qū)域v的磁化矢量；cα為第α個線圈敏感度矩陣；k為與梯度編碼磁場數目相同的向量，表示每個時刻下各個磁場凈梯度場；φ(x)為所有的梯度編碼磁場的多維度函數.

在非線性梯度場xy和線性梯度場x，y共同作用下，φ可以表示為

(4)

式中：a，b代表線性梯度場x，y強度，從式(4)可知，線性梯度場強度a，b決定了非線性梯度場中心的位置，其中常數項-ab可由靜磁場補償，即ΔB=ab，由此將上式進一步簡化為

(5)

施加線性梯度場后，非線性梯度場中心點從(0，0)轉移到位置(-a，-b)，當線性梯度場連續(xù)變化時，非線性梯度場中心則由一系列a，b值改變，經過多次移動，編碼磁場SEMs覆蓋整個成像空間，從而對整個成像平面編碼. 此方法與傳統的編碼策略不同之處在于，xy中心每移動至一個位置對成像空間編碼一次，近似于相位編碼操作，則k僅有讀出方向. 最后利用重建算法獲取重建圖像.

2 實驗方法

如上所述，由于編碼梯度場中心梯度低，非線性梯度成像中心模糊，故利用x，y梯度場用來移動非線性梯度場的中心，覆蓋整個成像區(qū)域.

圖2給出了xy中心隨線性梯度場強度a，b在圓形上變化的示意圖(中心由“o”標出). 為了更顯著觀察疊加SEMs中心變化，設定非線性梯度場幅度為1，線性梯度場強度a，b滿足a2+b2=0.5. 從圖中可以看到疊加的SEMs在FOV上變化，避免低梯度場區(qū)集中在某一成像區(qū)間.

設計非線性梯度場中心分布曲線，選擇圓形、雙曲線形、對角線形、螺線形. 為了在相同的條件下比較各種曲線的成像結果，仿真實驗中將各個SEMs強度定性的設為1(最大值). 式(5)中a，b分別滿足：

① 圓形：a2+b2=1；

② 雙曲線形：a2-b2=0.01；

③ 對角線形：a=±b；

④ 螺線形：a2+b2=l2，l∈[0，1].

確定最優(yōu)中心分布曲線形狀后，針對此形狀探討中心分布曲線大小，即線性梯度場強度變化對圖像重建的影響. 理論上高強度線性梯度場對圖像重建結果具有一定改善作用，但是同時會帶來PNS的風險，因此在尋找最優(yōu)成像效果的同時應盡量避免使用強度過高的線性梯度場.

最后，為驗證本文方法的有效性，在不同加速因子R下給出一組與O-Space技術對比實驗，加速因子是通過將中心移動次數減少到1/R次實現. O-Space技術的仿真實驗中，中心分布曲線根據文獻結果選擇圓形方案[10]. 所有的仿真都選擇8通道陣列接收線圈，線圈敏感度如圖3所示.

3 仿真實驗

仿真實驗和圖像重建在Matlab(The Mathworks， Natick， MA)中實現，計算機配置為32 bit，2.93 GHz雙核處理器，仿真使用Shepp-logan模型. 根據編碼過程，圖像重建實際上是解非對稱線性方程Ax=b，中心分布曲線仿真實驗中直接利用A的偽逆 (pinv(A))進行圖像重建. 在與O-Space對比試驗中，使用CG迭代算法重建[11]，根據圖像重建結果選擇合適迭代次數. 在偽逆重建中選擇了512×512編碼矩陣，迭代算法計算精度高，故選擇了128×128矩陣，在存在加速因子R且不為1的重建實驗中，矩陣尺寸由R值決定.

4 實驗結果

4.1 中心分布曲線

不同的中心分布曲線下xy梯度場圖像重建結果如圖4所示，為了便于比較，同時給出了Patloc技術重建圖像.

與參考圖像比較可以看到Patloc中心部分存在明顯模糊不清，在本文成像方法下，中心部分成像結果都有一定程度的改善，表明移動中心的方法能夠有效解決非線性梯度場中心成像模糊問題；對比各軌跡結果，圓形中心分布曲線重建結果最優(yōu)，較好地改善了中心模糊問題；雙曲線形中心分布曲線外側成像結果好，中心分辨率差，表明雙曲線形曲線在中心部分編碼不足；對角線形中心分布曲線結果中存在嚴重90°相位差偽影，xy非線性梯度場對稱軸沿對角線方向，故非線性梯度場中心沿對稱軸變化，90°相位差信號累加，造成圖像偽影；螺線型曲線中心點在FOV中心部分分布更為集中，重建結果中心結果差，分辨率明顯不如圓形軌跡結果. 對比結果表明xy非線性梯度場在圓形變化軌跡下對整個成像空間編碼最好.

設置圓形中心分布曲線a2+b2=r2，其中r為圓半徑，r從0開始增加. 重建結果進行歸一化，然后計算與參考圖像各像素點均方差(mean square error， MSE)之和. 根據MSE和峰值信噪比之間(peak signal to noise ratio，PSNR)的關系RPSNR=10lg(θmax2/θMSE)可知，MSE越高，對應的PSNR越低，圖像結果越差.

從圖5可知，r在[0， 1.5]段時，隨著線性梯度場增強，MSE顯著下降；當r=1.5時，MSE達到最低點，表示在此線性梯度場強度下非線性梯度場對成像空間編碼效率高，成像效果最好；進一步增加r，MSE上升，并在r>2.5之后，MSE保持在0.026 5附近，小范圍波動.

從上述分析可以得出結論：無限增加線性梯度場的強度對改善重建圖像無益，并且存在引起PNS的危險，非線性梯度場成像中合適的中心分布曲線應該根據實際需求選擇.

4.2 本文方法與O-Space對比實驗

為了提高圖像信噪比，利用CG迭代算法進行圖像重建，加速因子R分別選擇1，4，8.xy和Z2非線性梯度場幅度保持一致，Z2磁場中選擇z為0的平面. 根據上文結果，本文方法也選擇圓形軌跡. 在仿真數據重建中，當迭代次數超過8時，CG迭代算法收斂，所以在所有的仿真中都選擇8次迭代.

與直接求編碼矩陣的偽逆重建結果(圖4)對比顯示，迭代重建算法結果信噪比高很多，但時間消耗也較長，在上述計算機配置下完成8次迭代計算需要30 min. 與參考圖像對比表明，隨著加速因子R的增加，兩種方法重建結果都有所下降，但R=4時圖像結果MSE仍然較低，表明本文方法可以利用在快速磁共振成像中，提高成像速率；當R=8時，信噪比下降較快. 與O-Space技術相比，本文方法在所有的加速因子R下重建圖像MSE都更低，成像結果更好，表明xy非線性梯度場成像方法的有效性.

5 結 論

本文中研究了xy非線性梯度場成像技術. 通過仿真實驗，討論了xy非線性梯度場在線性梯度場作用下不同中心分布曲線的圖像重建結果，包括雙曲線、螺線形、圓形及對角線形曲線，結果表明圓形軌跡能夠得到最好的重建結果；然后針對圓形中心分布曲線，討論了線性梯度場強度，即圓半徑對成像的影響，在半徑r為1.5時得到重建結果最優(yōu)；迭代算法與直接偽逆相比表現出了優(yōu)越的重建結果，與O-Space的對比實驗中在所有的加速因子R下，本文方法圖像重建結果與參考圖像MSE更低，信噪比更高. 另外，Patloc研究組已制作雙曲線型非線性梯度場線圈[12-13]，能夠利用已有設備進行此方向的具體實驗驗證，下一步將針對xy非線性梯度場展開更多非線性梯度場并行成像方法，重建算法的研究.

[1] Hennig J， Welz A M， Schultz G， et al. Parallel imaging in non-bijective， curvilinear magnetic field gradients: a concept study[J]. MAGMA， 2008，21(1-2)：5-14.

[2] 何文勝.磁共振影像系統雙梯度技術的發(fā)展和現狀[J].中國醫(yī)療設備，2010，25:51-53.

He Wensheng. Development and current status of dual gradient technology in magnetic resonance imaging system[J]. China Medical Equipment， 2010,25:51-53. (in Chinese)

[3] Gallichan G， Cocosco C A， Dewdney A， et al. Simultaneously driven linear and nonlinear spatial encoding fields in mri[J]. Magn Reson Med， 2011,65(3):702-714.

[4] Stockmann J P， Galiana G， Tam L， et al. In vivo o-space imaging with a dedicated 12 cm z2 insert coil on a human 3t scanner using phase map calibration[J] Magn Reson Med， 2013,69(2):444-455.

[5] Stockmann J P， Ciris P A， Galiana G， et al. O-space imaging: highly efficient parallel imaging using second-order nonlinear fields as encoding gradients with no phase encoding[J]. Magn Reson Med， 2010,64(2):447-456.

[6] Tam L K， Stockmann J P， Galiana G， et al. Null space imaging: nonlinear magnetic encoding fields designed complementary to receiver coil sensitivities for improved acceleration in parallel imaging[J]. Magn Reson Med， 2012,68(4):1166-1175.

[7] 董薛.永磁型磁共振系統中勻場線圈設計的研究[D].武漢：華中師范大學，2011：11-13.

Dong Xue. Shim coils design in permanent magnet MRI system[D]. Wuhan: Huazhong Normal University， 2011:11-13. (in Chinese)

[8] Poole M. Improved equipment and techniques for dynamic shimming in high field MRI[D]. Nottingham: University of Nottingham, 2007.

[9] Schultz G， Ullmann P， Lehr H， et al. Reconstruction of mri data encoded with arbitrarily shaped， curvilinear， nonbijective magnetic fields[J]. Magn Reson Med， 2010,64(5):1390-1403.

[10] Stockmann J P， Ciris P A， Constable R T. Efficient o-space parallel imaging with higher-order encoding gradients and no phase encoding[C]∥Proc Intl Soc Mag Reson Med Honolulu: [s.n.], 2009:762.

[11] Pruessmann K P， Weiger M， Rnert P B， et al. Advances in sensitivity encoding with arbitrary k-space trajectories[J]. Magnetic Resonance in Medicine， 2001,46:638-651.

[12] Welz A M， Zaitsev M， Jia F， et al. Development of a non-shielded patloc gradient insert for human head imaging[C]∥Proc Intl Soc Mag Reson Med. Honolulu: [s.n.], 2009:3073.

[13] Cocosco C A， Gallichan D， Dewdney A J， et al. First in-vivo results with a patloc gradient insert coil for human head imaging[C]∥Proc Intl Soc Mag Reson. Montreal， Canada：[s.n.]， 2011:714.

(責任編輯：劉芳)

MRI Imaging Technology Based on xy Nonlinear Gradient

XU Ya-jie1，2， YANG Xiao-dong1

(1.Suzhou Institute of Biomedical Engineering and Technology， Chinese Academy of Sciences， Suzhou，Jiangsu 215163，China; 2.University of Chinese Academy of Sciences， Beijing 100049， China)

To eliminate center blurring existing in nonlinear gradient imaging, an imaging technology based onxynonlinear gradient field was proposed, and the effect of center displacement curve in magnetic resonance imaging (MRI) was discussed. The results reveal that circle center displacement curve with diameterr=1.5 achieves the best reconstruction, and in comparison experiment with O-Space under different acceleration factorR,xynonlinear magnetic field achieves a better image result.

nonlinear gradient imaging; center placement curve; CG iterative algorithm

2013-07-24

國家自然科學基金資助項目(11105096)；蘇州市科技資助項目(SYG201125，SH201207)

徐雅潔(1987—)，女，助理研究員，E-mail:xuyj@sibet.ac.cn.

楊曉冬(1977—)，男，研究員，E-mail:xiaodong.yang@sibet.ac.cn.

TP 391.41

A

1001-0645(2016)01-0064-06

10.15918/j.tbit1001-0645.2016.01.012