周蕾蕾，胡穎，郁蕓△，胡新華，楊坤

(1.南京醫(yī)科大學(xué)基礎(chǔ)醫(yī)學(xué)院，江蘇 南京 210029；2.南京醫(yī)科大學(xué)腦科醫(yī)院，江蘇 南京 210029)

1 引 言

腦功能磁共振成像(fMRI)是一種研究腦功能的非介入技術(shù)，且在空間上和時間上都有著很高的分辨率，已成為被廣泛使用的腦功能研究手段之一[1-2]。但大多數(shù)腦功能磁共振圖像均存在著噪聲污染[3]，給圖像分析以及臨床應(yīng)用造成了很大的困擾。為了適應(yīng)圖像數(shù)據(jù)的應(yīng)用需要，fMRI數(shù)據(jù)必須經(jīng)過一系列的預(yù)處理。磁共振圖像的預(yù)處理步驟在統(tǒng)計參數(shù)圖(SPM)中主要由以下幾個步驟構(gòu)成：片定時、重新排列、使標準化、平滑[4]。其中空間平滑技術(shù)的主要功能是消除被試在做磁共振時由于硬件不穩(wěn)定或被試生理運動產(chǎn)生的干擾信號，突顯圖像中的寬大區(qū)域和低頻成分、提高信噪比[5]。

目前國際上應(yīng)用于腦功能磁共振成像處理的軟件包有AFNI、SPM、FSL、MRIcro等。其中SPM技術(shù)已被fMRI和PET腦成像分析等領(lǐng)域廣泛采用，是腦功能研究不可或缺的技術(shù)。它是由英國倫敦大學(xué)的Friston教授等人基于Matlab平臺研究開發(fā)的，應(yīng)用隨機高斯場統(tǒng)計的方法來進行腦功能磁共振成像(fMRI)和PET的圖像數(shù)據(jù)的處理與分析的一個通用軟件包[6]。本研究基于SPM軟件包分析了圖像空間平滑技術(shù)的原理及其在處理fMRI數(shù)據(jù)的程序?qū)崿F(xiàn),并對其不足進行了改進。

2 SPM圖像空間平滑的算法原理及程序?qū)崿F(xiàn)

SPM圖像空間平滑技術(shù)采用具有低通濾波特性的高斯濾波器來進行平滑化處理，消除圖像重建過程中所產(chǎn)生的誤差以及不同被試的細微腦結(jié)構(gòu)差異所引起的誤差，達到提高圖像質(zhì)量的目的[7]。其基本原理是先進行圖像數(shù)據(jù)格式的轉(zhuǎn)換并確定參數(shù)最大半峰全寬(full width at half maximum，fwhm)，然后構(gòu)造卷積核函數(shù)，最后利用腦功能圖像與核函數(shù)進行卷積來實現(xiàn)圖像平滑。

2.1 卷積平滑

設(shè)P(x,y,z)為經(jīng)過校正、融合、分割、標準化后的圖像數(shù)據(jù)，G(x,y,z)為卷積核函數(shù)，Q(x,y,z)為經(jīng)過空間平滑處理后的圖像數(shù)據(jù)。

Q(x,y,z)=P(x,y,z)?G(x,y,z)

(1)

P(x,y,z)為經(jīng)過校正、融合、分割、標準化后的圖像數(shù)據(jù)。由于磁共振掃描儀記錄的功能圖像數(shù)據(jù)為DI-COM格式，而SPM中圖像處理的數(shù)據(jù)為IMG格式。所以，在運用SPM進行相關(guān)圖像處理前，需要進行圖像數(shù)據(jù)的格式轉(zhuǎn)換，可以利用相關(guān)軟件轉(zhuǎn)換工具來實現(xiàn)，如MRIConvert等。當(dāng)然，現(xiàn)階段的SPM的版本已經(jīng)自帶了一些數(shù)據(jù)格式轉(zhuǎn)換的功能。

(1)式用符號?來表示卷積運算，其數(shù)學(xué)意義是在定義區(qū)間內(nèi)，一個函數(shù)翻轉(zhuǎn)并平移后與另一個函數(shù)的乘積的積分。對于任意一個可積函數(shù)都可以構(gòu)造出一個逼近與其本身的光滑函數(shù)，且卷積得到的函數(shù)相比于原函數(shù)更加光滑。

某一方向上卷積的部分程序代碼實現(xiàn)：

if lx > 1;

for i=1:ly

u=X(:,i);

u=[flipud(u(1:Ex)); u; flipud(u([1:Ex]+lx - Ex))]；

U=sparse(conv(full(u),kx));

X(:,i)=U([1:lx]+2*Ex)

end

End

2.2 卷積核的構(gòu)造

在SPM中選擇了高斯函數(shù)作為卷積核函數(shù)。高斯函數(shù)的傅立葉變換等于高斯函數(shù)本身，即高斯函數(shù)的傅立葉變換頻譜是單瓣的。這一性質(zhì)使得高斯函數(shù)平滑的圖像可以很好地濾除不希望的噪聲信號并保留大部分有用信號。

以x軸方向的核函數(shù)的確定為例：

參數(shù)fwhm是指高度為最大處高度一半時譜帶的全寬。它很大程度上影響了平滑效果，其值越大，平滑程度就越好。在SPM中fwhm默認值為[8 8 8],在進行fMRI圖像空間平滑處理時，一般取值在4到10 mm之間[8]。

x軸方向的一維矩陣：x=[-s，-(s-1)，…，0，…，(s-1)，s]

其卷積核函數(shù)krn(x)的確定是利用B樣條插值并借助誤差函數(shù)來構(gòu)造的：

(2)

觀察式(2)可知，高斯濾波器是一個有權(quán)值的線性濾波器，每一點的像素值均被該點鄰域內(nèi)各點的像素加權(quán)均值所取替，且距離中心點的距離不同，相應(yīng)的權(quán)值也不同。

G(x)=Krn(x)/sum(Krn(x))

(3)

同理可以得到Y(jié)方向和Z方向上的一維矩陣：

y=[-s，-(s-1)，…，0，…，(s-1)，s]

z=[-s，-(s-1)，…，0，…，(s-1)，s]

同理可以得到Y(jié)方向和Z方向上的卷積核函數(shù)：

高斯函數(shù)是單值函數(shù)，這一特性決定了每一個鄰域像素點的權(quán)值隨該點與中心點的距離單調(diào)遞減。在平滑運算中與算子中心距離越大的像素點對算子中心的作用越小，這與圖像內(nèi)各像素點之間的作用關(guān)系是相吻合的，有利于圖像保真。

三是環(huán)境保護社會組織出現(xiàn)兩極分化的現(xiàn)象，一些有國際資金支持的組織、一些有國內(nèi)財團支持的社會組織和國際環(huán)境保護社會組織在華分支機構(gòu)在環(huán)境保護事業(yè)中迅速發(fā)展壯大，但一些本土的社會組織發(fā)展起色不大；環(huán)境保護社會組織數(shù)量仍然偏少，影響力總體仍然偏弱，建設(shè)性總體不足，作為全社會環(huán)境保護的參與和協(xié)調(diào)組織，難以填補政府、公民、中介技術(shù)服務(wù)組織和企業(yè)之間的角色空白，亟須立法予以經(jīng)濟、技術(shù)等方面的支持。

程序代碼實現(xiàn)：

if length(s)==1; s=[s s s]; end/*默認s是1×3的矩陣*/if isstruct(P),

VOX=sqrt(sum(P.mat(1:3,1:3).^2));/*調(diào)整系數(shù)的確定*/else

VOX=[1 1 1];

end;

s=s./VOX;

2.3 圖像數(shù)據(jù)與光滑函數(shù)卷積

可分離性是高斯函數(shù)的又一個重要性質(zhì)。鑒于這一特性，高斯函數(shù)空間卷積可以分兩步來做。首先將原圖像數(shù)據(jù)與一維高斯核函數(shù)進行卷積，然后將其卷積結(jié)果再與其方向垂直的一維高斯核函數(shù)進行卷積。

以x軸方向的圖像數(shù)據(jù)與光滑函數(shù)卷積為例，見圖1。

圖1卷積平滑示意圖

Fig1SchematicdiagramofConvolutionsmoothing

卷積平滑后圖像在x3處的像素值為：

x3′=(x1×1+x2×2+x3×4+x4×2+x5×1)/(1+2+4+2+1)

根據(jù)高斯函數(shù)的可分離性得：

Q(x,y,z) =P(x,y,z)?G(x,y,z)

=P(x,y,z)?G(x)?g(y)?g(z)

這樣，即實現(xiàn)了圖像空間平滑。

3 SPM平滑處理結(jié)果分析與改進

3.1 實驗結(jié)果分析

下列圖像為正常人腦的功能核磁共振成像的第46層的三視圖,且同組的1號圖為原始圖像，同組2號圖為經(jīng)過SPM空間平滑后的圖像。

圖2 SPM平滑前后的結(jié)果比較

Fig2Comparisonoforiginalimagesandsmoothedimages

對比圖2中a、b、c各組中1號圖像和2號圖像，不難看出經(jīng)過SPM圖像空間平滑后，圖像中的噪聲得到了很好的抑制，大部分有效信號得以保留，但也可以留意到圖像的部分邊緣存在模糊現(xiàn)象。結(jié)合文獻[9]和實驗結(jié)果，分析其原因可能是SPM圖像空間平滑時，部分有效低頻分量一定程度上被非線性增強，使得處理后的圖像存在部分邊緣模糊。

3.2 對SPM算法的改進及仿真

在一幅圖像中，有效信息對應(yīng)著低頻分量,噪聲對應(yīng)著高頻分量[10-11]。針對上述特點，本研究提出一種先進行一次預(yù)濾波，濾除一些頻率甚高的高頻分量，然后再進行SPM圖像空間平滑的改進算法,實現(xiàn)頻域濾波與空域濾波相結(jié)合的高效去噪且降低邊緣模糊的圖像平滑。

基本步驟：

(1)通過離散傅里葉變換將圖像由空域變換到的頻域(FFT)；

圖像中頻率表征了灰度變化劇烈程度，在頻率發(fā)生明顯突變處，多數(shù)有噪聲存在。

(2)預(yù)先濾除一些頻率甚高的高頻分量；

為了防止有效數(shù)據(jù)被錯誤濾除，只預(yù)先濾除一些頻率甚高的高頻分量，即只濾除了部分高頻噪聲分量，故仍有相當(dāng)一部分次高噪聲存在。

(3)離散傅里葉逆變換(IFFT)；

將經(jīng)過(2)預(yù)處理后的圖像數(shù)據(jù)變換到空域。

(4)SPM空間平滑。

沿用SPM中以高斯核函數(shù)為模版的空間卷積平滑，實現(xiàn)噪聲圖像的濾波去噪處理。

改進后算法流程見圖3。

圖3改進后算法的平滑框圖

Fig3Blockdiagramoftheimprovedalgorithm

預(yù)濾波過程：

為了能更好的減小圖像的模糊程度，更多的保留圖像邊緣信息。且考慮到在步驟(4)中會再一次進行濾波處理。故本研究中選取截止頻率相對較大的理想低通濾波器來進行圖像的預(yù)濾波處理。

R(u,v)為圖像中的某一點(u,v)到中心點的距離，R0為截止頻率。

本研究以3.1組中的c1原始圖像為例進行實驗，利用改進后的SPM圖像空間平滑算法進行圖像平滑化處理，得到如圖4中所示的效果對比圖 。

圖4 不同方法平滑后的結(jié)果比較

3.3 實驗結(jié)果分析

圖5不同方法平滑后的圖像灰度直方圖的結(jié)果比較

Fig5Comparisonofgrayhistogramofdifferentsmoothedimages

圖5給出的是與圖4中的原始圖像、SPM平滑后的圖像以及改進后的平滑圖像相對應(yīng)的灰度分布直方圖。

結(jié)合圖4中不同方法平滑后的效果圖和圖5中的不同方法平滑后的圖像灰度分布直方圖，比較改進后的SPM圖像空間平滑結(jié)果與原SPM的圖像空間平滑結(jié)果，發(fā)現(xiàn)改進后的SPM圖像空間平滑方法不僅很好的降低了圖像中的噪聲，也減少了圖像的邊緣模糊。使原SPM圖像空間平滑中存在的不足得到了改善。

4 總結(jié)

腦功能磁共振成像技術(shù)探測到的神經(jīng)生理信號來源于數(shù)毫米尺度內(nèi)腦血流等的改變,圖像重建中有效信號對應(yīng)著低頻率部分,而噪聲對應(yīng)著高頻部分[12]。SPM采用一種具有低通濾波特性的高斯濾波器進行平滑，利用高斯分布函數(shù)構(gòu)造高斯模板圖像,將腦功能圖像和模板圖像進行卷積,來實現(xiàn)了平滑化濾波處理。

由3.1中各組2號實驗結(jié)果顯示可知，SPM圖像空間平滑后噪聲得到很大抑制，減少了圖像重建過程中所產(chǎn)生的誤差和不同被試腦結(jié)構(gòu)之間的細微差別。但其不足之處也比較明顯，在高效濾去高頻噪聲分量的同時對低頻有效信號分量存在著一定程度的非線性增強，使處理后的圖像因?qū)Ρ榷炔幻黠@而存在模糊。本研究通過改進SPM平滑算法，一定程度上減少了圖像模糊，提高了圖像的質(zhì)量。

總之，大腦是人體的一個重要器官，也是一個結(jié)構(gòu)比較復(fù)雜的組織，腦部圖像處理是醫(yī)學(xué)圖像處理的重要研究方面，對臨床診斷與治療有著重要意義。腦功能磁共振成像技術(shù)在實際采集圖像過程中會產(chǎn)生很多不同種類不同性質(zhì)的噪聲，SPM圖像空間平滑是圖像處理中相對較好的圖像處理技術(shù)。為了適應(yīng)更精準更科學(xué)的醫(yī)療診斷需求，它仍然存在著很大的改進空間。