融合NSST和稀疏表示的PET和MRI圖像融合

康家銀，陸 武，張文娟

1(江蘇海洋大學 電子工程學院，江蘇 連云港 222005) 2(連云港市第一人民醫院 核醫學科，江蘇 連云港 222061) 3(江蘇海洋大學 計算機工程學院，江蘇 連云港 222005)

1 引 言

正電子發射斷層(Positron Emission Tomography，PET)是一種核醫學成像技術，也是當今最先進的無創傷性診斷手段之一[1].磁共振成像(Magnetic Resonance Imaging，MRI)是一種能夠呈現人體解剖信息的結構型成像[2].PET圖像具有能夠反映人體不同組織的新陳代謝情況的顏色信息，但其空間分辨率低，結構信息不明顯；MRI圖像的空間分辨率高、軟組織結構信息強，但缺乏反應人體新陳代謝功能的顏色信息[3,4].通過影像融合技術可以將PET圖像和MRI圖像各自的優勢有效結合，從而使得融合后的圖像在保留PET圖像顏色信息的同時，又具有MRI圖像明顯的空間結構信息，以便融合后的圖像將更有助于臨床中對疾病的診斷和治療[5].

近年來，針對PET和MRI圖像的融合，研究者們已經提出了不同的策略和算法.其中一類是基于變換和替換策略的IHS(Intensity-Hue-Saturation)方法[5,6].這類方法首先將PET變換到IHS空間中；然后用匹配后的MRI(將MRI與PET在IHS空間中的I分量進行直方圖匹配)去替換PET在IHS空間中的I分量；最后經反變換得到融合后的PET圖像.這類方法通常會得到一幅空間信息強(分辨率高)的PET圖像，但會使得在信息替代的過程中損失PET圖像的顏色信息.另一類用于PET和MRI圖像融合的方法是基于多分辨率分析(Multi-Resolution Analysis,MRA)的方法，如小波變換等[7,8].這類方法首先將PET和MRI圖像進多尺度分解(分解成系數和變換基)；然后將分解后的系數根據某種融合規則進行融合；最后將融合后的系數和變換基進行反變換，從而得到最終的融合圖像.這類方法可以有效的保留PET圖像的色彩信息，但在提升PET圖像的空間結構信息方面具有不可避免的局限性.第三類用于PET和MRI圖像融合的方法是基于稀疏表示(Sparse Representation，SR)的方法[9,10].這類方法首先求出PET和MRI圖像各自的稀疏表示系數；然后將二者的稀疏系數采用某種融合規則進行融合；再利用融合后的稀疏系數和用于重構的過完備字典重構出最終的融合圖像.稀疏表示在圖像融合上取得了引人注目的效果，但目前用于PET和MRI圖像融合的稀疏表示模型中,用于圖像重構的字典的訓練/構造是針對整幅圖像進行的，即從整幅圖像中提取圖像塊來構造全局字典，從而導致計算復雜度的增加.此外，由于構造全局字典時沒有考慮圖像中塊與塊之間的結構相似性，致使求解得到的稀疏系數不夠準確，從而難以準確的重建出目標圖像[11].

非下采樣剪切波變換(Non-Subsampled Shearlet Transform，NSST)[12]，一種先進的多分辨率、多尺度分析方法.除了具有其它一些多分辨率分析方法具有的局部性、多尺度和多方向性特點外，NSST還具有平移不變性、計算復雜性低等特點.因此，NSST被成功的應用于圖像融合中[13-15].

根據以上分析，本文結合NSST和稀疏表示各自的優點，提出一種基于NSST和稀疏表示的PET和MRI圖像融合方法和框架.首先，利用NSST將PET和MRI圖像分解成高頻和低頻成分.然后，使用梯度值較大法將對應的高頻部分進行融合；基于SR將對應的低頻部分進行融合，即利用SR獲取能表征低頻成分的稀疏系數，再對獲取的稀疏系數進行融合，并獲得融合的低頻成分.最后，利用NSST逆變換得到PET與MRI的融合圖像.該方法不僅能夠有效的保持PET圖像中的顏色信息，而且還能夠很好的集成MRI圖像中豐富的結構信息.

2 相關工作

2.1非下采樣剪切波變換

不同于之前提出的多分辨率分析(MRA)方法(如輪廓波變換[16]、曲線波變換[17]等)，非下采樣剪切波變換(NSST)是一種新穎的多分辨率、多尺度分析方法.NSST除具有其它MRA方法的特點，如局部性、多尺度和各向異性等外，還具有計算復雜度更低、實現效率更高、分解方向數不受限制等優點[13].

圖1 NSST分解示意圖Fig.1 Diagram of NSST decomposition

2.2 圖像的稀疏表示理論

近年來，稀疏表示(SR)成為研究的熱點并已成功的應用到圖像處理和計算機視覺的各個方面，如超分辨率重建[18]，圖像重構[19]，多源/多傳感器/多聚焦圖像的融合[20-22]等.圖像的稀疏表示模型是基于一種假設，即一幅圖像可以被事先定義好的過完備字典中的“原子”稀疏性的表示[23].稀疏表示使得圖像的信息主要集中在少數“原子”上，而這些“原子”反映了圖像的主要特征和內部結構.圖像的稀疏表示模型可表示為：

(1)

其中x∈n是大小為的圖像塊(以列向量形式表示)；α∈N是稀疏表示系數；‖·‖0是l0范數，表示向量α中非零元素的個數；D∈n×N(n?N)是過完備字典；ε是預先設定的逼近誤差容限，且ε≥0.

式(1)是一個非凸函數，所以上述稀疏表示最優化問題，即(1)的求解是一個NP-難問題[24].因此需要采用特殊的方法對(1)進行求解.其中一類常用的方法是用l1范數代替l0范數，從而將(1)中的非凸優化問題轉化為l1范數的凸優化問題，于是便可得到如下的基追蹤(Basis Pursuit，BP)優化問題：

(2)

3 本文算法

3.1 算法框架

基于前面引入的NSST和SR基本原理，本文提出的PET和MRI圖像融合基本框架如圖2所示(假設源PET和MRI圖像已經進行了配準)，其實現步驟為：

2)將PET和MRI高頻子帶采用梯度值最大規則進行融合，具體的融合規則將在隨后的小節中詳述.

3)對于二者的低頻子帶利用SR進行融合.類似的，融合規則的具體實現方法將在后面的小節中給出.

4)最后，將整合后的高、低頻子帶進行NSST反變換，從而得到最終的融合圖像.

圖2 PET和MRI圖像融合框架Fig.2 Image fusion framework for PET and MRI images

3.2 基于梯度值法的高頻子帶融合

圖像經過NSST分解后，邊緣、輪廓等許多源圖像的細節信息包含在高頻子帶中.高頻子帶融合的目的就是從源圖像中盡可能多的保留這些有用的細節信息.傳統的基于“絕對值取大”的高頻子帶融合規則僅僅考慮了子帶圖像中“像素”自身，而沒有考慮該像素與其周圍像素間的空間關系.因此簡單的采用傳統的基于“絕對值取大”的高頻子帶融合規則容易導致源圖像中的細節信息的丟失和引入額外的人為 “噪聲”.為此，本文采用如下的基于“梯度值取大”的高頻子帶融合規則：

(3)

其中IPH(i,j)，IMH(i,j)和IFH(i,j)分別為PET，MRI和融合的高頻子帶圖像中(i,j)處的值；grad(·)為求PET和MRI高頻子帶中(i,j)處梯度值的操作，本文采用Sobel算子求梯度值，其形式如圖3所示.

圖3 Soble梯度算子Fig.3 Sobel gradient operator

3.3 基于SR的低頻成分融合

圖像經過NSST分解后，源圖像的信息主要集中在低頻子帶圖像中.因此，低頻子帶的融合非常關鍵.為此，本文基于稀疏表示進行PET和MRI低頻子帶圖像的融合.

3.3.1 自適應局部字典的構造

在圖像的稀疏表示模型中，字典(無論是用于稀疏系數求解的編碼字典還是用于圖像表示的重構字典)的構造至關重要.經典的稀疏表示中，字典的構造或學習都是基于整幅圖像的，即構造或學習得到的是針對整幅圖像的全局字典，從而難以體現圖像中的局部相似性.此外，全局字典的引入也勢必會增加計算復雜度.因此，本文采用如下的策略構造能夠反應圖像中塊與塊之間結構相似性的局部自適應字典.

圖4 局部字典構造示意圖Fig.4 Diagram of local dictionary construction

3.3.2 低頻子帶系數的融合

本文提出的基于稀疏表示的PET和MRI圖像低頻子帶融合框架如圖5所示，分為四個階段(見圖5中的①，②，③，④)，具體的實現過程如下：

第一階段：聯合字典的構建

第二階段：稀疏編碼

第三階段：稀疏系數的融合

圖5 基于SR的低頻子帶融合示意圖Fig.5 Diagram of low frequency component fusion based on SR

第四階段：稀疏重構

重復上述1)到7)步，直到將低頻子帶(IPL和IML)中所有的“像素點”(共N對)處理完畢.

4 實驗結果與分析

為了驗證本文所提出的用于PET和MRI圖像融合的算法框架的性能，我們用一組公開的數據集進行了相關的實驗.

4.1 實驗數據介紹

本文用于實驗的PET和MRI圖像來自于哈佛大學醫學院的“全腦”(The Whole Brain)數據集1.本文選用了30組PET和MRI圖像進行了融合實驗，圖像大小均為256×256.

4.2 算法參數設置

本文實驗中算法的參數設置如下：

NSST分解級數：L=3；

搜索窗口滑動步長：6；

“最”相關的圖像塊個數：s=60.

4.3實驗結果分析

本文選取了30組數據進行實驗，并將本文方法得到的結果與其它方法的結果進行了比較.其中用于比較的方法分別來自于文獻[5](記為IHS-Retina)、文獻[6](記為IHS-Gabor)和文獻[10](記為CT-SR).圖6和圖7為兩組數據(分別以Coronal和Transaxial兩個視角呈現)的融合結果.

圖6 PET和MRI圖像的融合結果(Coronal視角)

圖7 PET和MRI圖像的融合結果(Transaxial視角)

由圖6和圖7可知，就保留MRI圖像的結構信息而言，算法NSST(融合結果分別見圖6(e)和圖7(e))和本文算法(融合結果分別見圖6(h)和圖7(h))明顯好于其它幾種方法，見圖6和圖7的各個子圖中用矩形框 "標注"的區域；而在其它方法中，IHS-Retina方法的效果最不理想，即模糊化程度最嚴重(見圖6(c)中的橢圓形區域和圖7(c)中的圓形區域).此外，相較于圖6(e)和圖7(e)，圖6(h)和圖7(h)的結構信息更接近原始的MRI圖像，見圖6和圖7的各個子圖中用橢圓形框 "標注"的區域.

1http://www.med.harvard.edu/aanlib/

另外，由圖6和圖7可知，在保留PET圖像顏色信息方面，本文方法略遜于IHS-Retina方法，與其它方法相當.此外，將各融合結果與原始MRI圖像(圖6(b)、圖7(b))的腦殼(Skull)部分進行對照，可見唯有本文方法得到的融合結果與原始MRI圖像最接近(灰度一致)，而其它方法得到的結果中均有不同程度的偽彩色.

綜合以上兩方面的因素可知，就視覺效果而言，相對于其它幾種方法，本文方法取得了更好的融合效果，即在較好的保持PET圖像顏色信息的同時，很好的融合了MRI圖像的結構信息.

為了進一步評價融合效果，本文采用一些客觀評價指標對融合結果進行客觀、定量的評價和分析.本文用于評價融合效果的指標包括信息熵(Information Entropy，IE)[25]、平均梯度(Average Gradient，AG)[25]、光譜差異性(Spectral Discrepancy、SD)[5]和Qabf[26].各個評價指標的具體形式如下：

信息熵IE主要用于衡量融合后圖像中所含信息量的大小.IE越大，表明圖像中包含的信息越豐富.本文中，IE越大，表明融合后，PET圖像包含的信息(色彩信息和結構信息)越豐富.IE定義如下：

(4)

其中P(i)表示灰度值為i(i=0,1,…,L-1)的像素出現的概率.

平均梯度AG反應的是融合圖像的清晰度(或空間分辨率)的高低.AG越大，表明圖像層次越多，越清晰.本文中，AG越大，表明融合后，PET圖像層次越豐富，空間分辨率越高.AG定義如下：

k=R,G,B

(5)

光譜差異性SD衡量的是融合圖像的光譜(顏色)特性；具體的，SD反應了融合圖像與源PET圖像的顏色差異性.SD越小，表明融合圖像的光譜特性越好.本文中，SD越小，表明融合后，PET圖像的色彩信息和原始PET圖像越接近.SD定義如下：

(6)

其中Fk(x,y)和Ok(x,y)分別為融合后的圖像與輸入的PET圖像在(x,y)處的值；R，G，B為PET圖像的三個顏色分量.

Qabf主要用于評價融合過程中源圖像的邊緣信息保持程度.Qabf值越大，表明融合圖像與源圖像間的邊緣相似度越高，或源圖像的邊緣信息損失越小.本文中，Qabf值越大，主要表明融合后的PET圖像中包含了更多的MRI圖像的邊緣信息(結構信息).Qabf定義如下：

(7)

利用不同方法進行PET和MRI圖像融合后得到的量化指標如表1所示(其中表1中的結果是30組數據的平均值).

表1 不同方法融合PET和MRI圖像的結果
Table 1 Fusion results of PET and MRI images using different methods

融合方法IEAGSDQabfIHS-Retina5.02210.1175.9070.565IHS-Gabor5.59710.9136.4310.635NSST5.66311.0856.2740.647SR5.69711.5066.1370.612CT-SR5.74110.4027.5510.563本文方法5.83612.0516.0290.681

注：加粗者為最優結果.

由表1可知，就信息熵(IE)、平均梯度(AG)和Qabf而言，本文方法的性能均優于其它方法.就光譜差異性(SD)而言，本文方法的性能略遜色于IHS-Retina，但優于其它四種方法.綜上可知，就定量指標而言，相對于其它幾種方法，本文提出的方法取得了更好的融合效果.

4.4 算法參數L對融合性能的影響

本文提出的算法框架中涉及到一個重要的參數，即NSST的分解層級L.本小節結合4.3節中給出的四個客觀量化指標(IE、AG、SD和Qabf)分析不同的L對融合性能的影響.圖8中，(a)、(b)、(c)和(d)分別為不同分解層級L對應的IE、AG、SD和Qabf.

圖8 分解層級L對融合性能的影響

由圖8(a)和圖8(d)可知，隨著分解層級的加大，融合性能會隨之提升，但當L>3時，性能提升便不明顯.同樣的，由圖8(b)和圖8(c)可知，隨著分解層級的加大，融合性能會隨之提升，但是提升不明顯.此外，分解層級越大，運算耗時便越長.因此，綜合考慮上述各個因素，最終選擇了L=3的分解層級.

5 結 論

本文提出了一種用于PET和MRI圖像融合的方法，即利用非下采樣剪切波變換(NSST)和稀疏表示(SR)的PET和MRI圖像融合方法.該方法充分的結合了NSST和SR各自的優點.利用公開的腦部數據集進行實驗的結果表明，本文提出的方法取得了很好的融合效果，即在有效保留PET圖像“顏色”信息的同時，很好的保持了MRI圖像的空間結構信息.此外，與其它已有的用于PET和MRI圖像融合的方法相比，本文的方法無論在視覺效果還是在客觀的定量評價中，都取得了最好的圖像融合效果.