基于互信息準(zhǔn)則的圖像平滑和分割

溫鐵祥潘正洋辜 嘉(中國科學(xué)院深圳先進(jìn)技術(shù)研究院深圳58055)

2(中國科學(xué)院大學(xué)北京100049)

3(中國科技大學(xué)合肥230000)

基于互信息準(zhǔn)則的圖像平滑和分割

溫鐵祥1,2潘正洋1,3辜 嘉1
1(中國科學(xué)院深圳先進(jìn)技術(shù)研究院深圳518055)

2(中國科學(xué)院大學(xué)北京100049)

3(中國科技大學(xué)合肥230000)

在計(jì)算機(jī)視覺領(lǐng)域，尺度空間扮演著一個(gè)很重要的角色。多尺度圖像分析的基礎(chǔ)是自動(dòng)尺度選擇，但它的性能非常主觀和依賴于經(jīng)驗(yàn)。基于互信息的度量準(zhǔn)則，文章提出了一種自動(dòng)選取最優(yōu)尺度的模型。首先，研究專注于基于形態(tài)學(xué)算子的多尺度圖像平滑去噪方法，這種技術(shù)不需要噪聲方差的先驗(yàn)知識，可以有效地消除照度的變化。其次，通過遞歸修剪Huffman編碼樹，設(shè)計(jì)了一個(gè)基于聚類的無監(jiān)督圖像分割算法。一個(gè)特定的聚類數(shù)從信息理論的角度來看，提出的聚類算法可以保留最大的信息量。最后，用一系列的實(shí)驗(yàn)對算法的性能進(jìn)行了驗(yàn)證，并從數(shù)學(xué)上進(jìn)行了詳細(xì)的證明和分析，實(shí)驗(yàn)結(jié)果表明本文提出的算法能獲得最優(yōu)尺度的圖像平滑和分割性能。

尺度選擇；高帽變換；互信息；形態(tài)學(xué)

1 引 言

尺度空間在許多領(lǐng)域的數(shù)學(xué)建模過程中是一個(gè)很重要的概念，但尺度的選擇非常復(fù)雜。在各研究領(lǐng)域，如：化學(xué)、物理、生態(tài)系統(tǒng)等，宏觀的連續(xù)行為產(chǎn)生于微觀粒子(分子、細(xì)胞、人群中的個(gè)體)之間的相互作用以及它們與周圍環(huán)境的之間的相互作用。在許多情況下，可通過使用各種宏觀模型的推導(dǎo)(仿真、優(yōu)化、分岔分析)對一個(gè)宏觀模型(流體流動(dòng)的Navier-Stokes方程或反應(yīng)擴(kuò)散方程)進(jìn)行定量的描述[1–3]。然而，對于許多復(fù)雜的系統(tǒng)，盡管演化是在宏觀尺度變換下觀察的，但是準(zhǔn)確的模型只能在給出更精細(xì)尺度和微觀模型上進(jìn)行描述，如：格子玻爾茲曼、活性粒子動(dòng)力學(xué)理論、分子動(dòng)力學(xué)、元胞自動(dòng)機(jī)等[4–6]。最近，在物理化學(xué)建模方面，宏觀和微觀的聯(lián)系已取得了初步的成效[7–9]。

尺度空間理論的應(yīng)用非常廣泛，而尺度空間中尺度之間的聯(lián)系卻是復(fù)雜的。本研究只關(guān)注在圖像處理任務(wù)方面的最優(yōu)尺度選擇。尺度分析是許多圖像處理任務(wù)方面最有用的框架之一[10–12]。分層多尺度表達(dá)的基本思想是分析有關(guān)尺度變換[13,14]，而它通常被定義為被過濾后的圖像集合。

在一些經(jīng)典文獻(xiàn)中，多尺度表示可以區(qū)分為三種類型：第一類是小波變換[15]，它在圖像處理領(lǐng)域被廣泛應(yīng)用，如：圖像去噪[16]、圖像和視頻壓縮[17,18]及多模態(tài)圖像配準(zhǔn)[19]。通過采用可變大小的窗口技術(shù)，小波變換可以把一個(gè)信號分解成多尺度表示；第二類是基于偏微分方程的擴(kuò)散過程[20]，它源于傳統(tǒng)的各向同性熱流動(dòng)方程，并已廣泛使用在圖像濾波[21–24]，分割[25–27]和圖像修復(fù)[28–30]；第三類是基于數(shù)學(xué)形態(tài)學(xué)[31–33]。大多數(shù)形態(tài)過濾器產(chǎn)生非線性尺度，與線性尺度分析相比，它在絕大多數(shù)圖像處理任務(wù)方面已經(jīng)表現(xiàn)出了卓越的性能。如：一幅圖像的特征(邊緣、角落等)可以很好的用非線性算子保留起來。

尺度是多尺度分析應(yīng)用中的關(guān)鍵參數(shù)，它的性能好壞決定了多尺度算法的成功與否。最近，一些現(xiàn)有的選擇合適尺度去除圖像噪聲的技術(shù)已經(jīng)被提出，這些技術(shù)通過最優(yōu)控制理論[34,35]、奇異攝動(dòng)方法[36]、最小化相關(guān)信號和噪聲[37]、統(tǒng)計(jì)模型[38,39]和馬氏隨機(jī)場模型[40]獲得最優(yōu)的尺度。然而，這些方法需要預(yù)先知道噪聲方差，且它的評估運(yùn)算非常耗時(shí)，有時(shí)也不切實(shí)際。為了提高尺度選擇的魯棒性，本文基于信息理論模型[41]，提出了一種新的尺度選擇的算法，該算法是基于最大化互信息來表示形態(tài)學(xué)多尺度。

此外，本文提出的尺度選擇模型可以用于無監(jiān)督圖像分割，其中的聚類數(shù)(Number of Cluster, NC)的確定被認(rèn)為是最優(yōu)尺度選擇。測量的標(biāo)準(zhǔn)來源于互信息的計(jì)算?；バ畔?biāo)準(zhǔn)已經(jīng)在基于二進(jìn)制空間的圖像分割方法[42]、模擬退火法[43]和模糊C-均值法(FCM)[44,45]中被采用。不同于那些方法，本文所提出的算法可以保持與原始圖像和分割后的圖像之間的最大信息量，并且模型解的凸性也可以通過采用哈夫曼編碼策略加以保證，并在數(shù)學(xué)上對算法良好的性能進(jìn)行了證明。

本文第2部分詳細(xì)描述了建立在互信息最大化目標(biāo)函數(shù)的基礎(chǔ)上的尺度選擇模型和它在圖像去噪和分割方面的應(yīng)用。第3部分介紹實(shí)驗(yàn)結(jié)果，第4部分對結(jié)果進(jìn)行討論。

2 互信息準(zhǔn)則模型

2.1 模型的描述

在信息理論中，互信息被定義為一種統(tǒng)計(jì)學(xué)上對兩個(gè)隨機(jī)變量X和Y的相關(guān)性測量。將圖像u(x)看作一個(gè)隨機(jī)變量，則u(x)的概率密度可表示為pu(i)，并且通?？蓮膱D像直方圖中可以估計(jì)出來。那么，給定兩幅圖像u(x)和v(x)，它們間的互信息I(u,v)可以寫為：

其中，H(u)表示香農(nóng)熵，表示一個(gè)隨機(jī)變量的平均信息或不確定性，并且只有在同等概率的情況下，H(u)表示的信息量達(dá)到最大。H(u,v)表示的是u(x)和v(x)的聯(lián)合熵，表達(dá)了這兩個(gè)變量之間的匹配性。puv(i,j)是衡量u(x)和v(x)之間的相似性的聯(lián)合概率。聯(lián)合概率可以定義為：

在隨機(jī)變量u(x)和v(x)正確匹配的情況下，互信息量達(dá)到最大化。為了表示原始圖像u(x)和用多度算法處理后圖像之間的關(guān)系，我們定義一個(gè)基于信息理論的能量函數(shù)Jt如下：

其中，Pt是尺度為t的一個(gè)操作算子，如：用于圖像去噪的形態(tài)學(xué)開操作算子和用于圖像分割的聚集操作算子。

通過尋找多尺度 t，最大化互信息函數(shù)：

原始圖像和處理后的圖像之間的互信息可表示為：

原始圖像和處理后的殘余圖像之間的互信息計(jì)算可表示為：

隨著尺度變量t的增長，可以看出來公式(7和(8)之間的互信息量的增長方向是相反的：一個(gè)在增長的同時(shí)另一個(gè)在下降。為了便于比較公式(7)和(8)，我們歸一化并重寫為：

在本文下面的章節(jié)中，分別用上面化一歸的方程作為圖像去噪和分割上問題上尺度選擇的標(biāo)準(zhǔn)。

2.2 基于互信息最大化的圖像去噪

2.2.1 去噪算法的描述

這一節(jié)將介紹基于多尺度表示的形態(tài)學(xué)開操作算子的圖像去噪算法，用具體的開操作算子Ot取代公式(9)和(10)中的尺度算子。開操作是膨脹和腐蝕操作的結(jié)合，是數(shù)學(xué)形態(tài)學(xué)中的基礎(chǔ)操作。其中，結(jié)構(gòu)元素形狀和尺寸大小的選擇非?；A(chǔ)和重要。然而，對于一幅具體圖像的處理任務(wù)，結(jié)構(gòu)元素的選擇往往依靠主觀和經(jīng)驗(yàn)。為了自動(dòng)地選擇結(jié)構(gòu)元素，我們從多尺度的觀點(diǎn)出發(fā)，用尺度t參數(shù)化結(jié)構(gòu)元素的大小并且重新定義基礎(chǔ)的形態(tài)算子。這樣，對于一個(gè)具體的圖像處理任務(wù)，選擇結(jié)構(gòu)元素的尺寸問題便成了選擇合適大小或者最優(yōu)尺度的問題。

膨脹算子D計(jì)算的是在給定鄰域的灰度最大值，結(jié)構(gòu)元素用尺度t進(jìn)行參數(shù)化，膨脹算子可表述為：

當(dāng)圖像被膨脹處理后，那些比結(jié)構(gòu)元素更小的圖像細(xì)節(jié)(如：噪聲等)將從圖像中消除。當(dāng)尺度t從0開始增大時(shí)，能獲得一系列過濾后的圖像。

相反地，腐蝕因子E計(jì)算鄰域附近的灰度最小值，定義為：

然后開操作算上O可以被定義為在腐蝕運(yùn)算后再進(jìn)行膨脹運(yùn)算：

用泰勒展式展開公式(13)，可得：

其中，o(t)是關(guān)于尺度t的高階無窮小量，▽u(x)是原始圖像u(x)的梯度。因此，隨尺度t的減小，開操作運(yùn)算的極限為：

同時(shí)，隨尺度t的增大，開操作運(yùn)算的極限為：

其中，C是常數(shù)。公式(16)表示：如果選定的尺度足夠大，通過開操作運(yùn)算得到的濾波后圖像是一個(gè)常量。

原始圖像與開操作運(yùn)算的殘差恰好為高帽變換：

開操作算子可以有效地抑制圖像中的細(xì)節(jié)信息，因此，當(dāng)開操作圖像從原始圖像減去后，就可以獲得所需的細(xì)節(jié)[46]。高帽變換的其中一個(gè)潛在應(yīng)用是估量不均勻的背景光照。然而，如果尺度過小，高帽變化會(huì)對不均勻光照估計(jì)不足。反之，高帽變化會(huì)對不均勻光照估計(jì)過度。也即，如果選擇的尺度太小的話，從高帽變換得到的結(jié)果會(huì)對噪聲很敏感。另一方面，如果選擇的尺度太大的話，高帽變換對原始的圖像起不了作用，因?yàn)閷Υ蟪叨葧r(shí)開操作的結(jié)果是一個(gè)常量圖像。利用公式(16)和(17)，以上結(jié)論可用數(shù)學(xué)公式描述為：

為了給基于高帽變換的圖像去噪選擇合適的尺度，利用互信息公式(9)和(10)作為定量標(biāo)準(zhǔn)，給出如下的尺度選擇算法(i-De)。

Algorithm:i-De

Parameter:

Denote T as the number of iteration, at the initial scale t＝1.

Do fort＝2, …, T

1.Compute top-hat transformation as in (17);

2.Calculate mutual information I1and I2as follows: as a small threshold parameter.

Calculate

2.2.2 i-De算法的分析

對于i-De算法的可行性，將對其中公式的單調(diào)性給出相關(guān)證明。下面首先介紹與開運(yùn)算相關(guān)的定理。

定理1：令為一函數(shù)集，并且它在對比度變化的情況下保持穩(wěn)定，另外O是一個(gè)具有sup-inf運(yùn)算形式的算子，定義為：那么，sup-inf算子是單調(diào)的并且在對比不變的。

根據(jù)定理1給出的開操作算子的單調(diào)性，下面介紹兩個(gè)描述公式(20)和(21)性質(zhì)的命題。

命題1：令O是一個(gè)開操作算子，T是一個(gè)高帽變換，I(u,v)是兩個(gè)隨機(jī)變量u(x)和v(x)的互信息。那么，公式(20)定義的互信息是單調(diào)遞減的，公式(21)定義的互信息是單調(diào)遞增的。

證明：公式(15)和(16)表明開操作后的圖像在尺度為0時(shí)是原始圖像，在尺度為無窮大時(shí)趨于常量圖像。結(jié)合定理1中的單調(diào)性，可以得出：隨著尺度從零增加到無窮大，開操作后圖像Otu(x)將從原始圖像單調(diào)下降為常量圖像。同時(shí)，因?yàn)樽畲蠡バ畔⒅挥挟?dāng)圖像與自身進(jìn)行互信息計(jì)算時(shí)才能獲得，即。所以，可以得出結(jié)論：公式(20)的互信息是單調(diào)減少的。公式(21)的互信息是單調(diào)增加的。

在實(shí)驗(yàn)結(jié)果中，觀察到公式(20)的互信息為單調(diào)減少，這和命題1的結(jié)論相吻合。然而，公式(21)的互信息為微小振蕩地緩慢增加，但不是單調(diào)的。這一現(xiàn)象是由于公式(20)的計(jì)算對噪聲不敏感。因?yàn)殚_操作圖像可以看作是原始圖像的一個(gè)平滑去噪過程。相反地，高帽變換圖像是原始圖像的開操作圖像的殘留圖像。也即，高帽變換圖像在小尺度時(shí)對于噪聲敏感，正如公式(18)所示。這樣的振蕩是由噪聲的存在所造成。但是公式(21)的互信息的趨勢是增加的。

命題 2：公式(20)和(21)定義的互信息將會(huì)在某個(gè)尺度處相交。

證明：通過將公式(15)和(16)替換成(20)，我們可以得出結(jié)論：公式(20)的互信息從最大的互信息(即1)減少到最小的互信息(即0)。相似地，通過將公式(18)和(19)替換成(21)，我們可以得出結(jié)論：公式(21)從0增加到1。因此，公式(20)和(21)會(huì)在某個(gè)尺度處相交。

根據(jù)命題2，隨著尺度t的增加，曲線(20)和(21)有一個(gè)相交點(diǎn)，并且這個(gè)交點(diǎn)是一個(gè)拐點(diǎn)，因?yàn)檫@兩個(gè)互信息在此開始收斂。然而，正如在我們實(shí)驗(yàn)結(jié)果中表明的那樣，該點(diǎn)的尺度不是最優(yōu)的。i-De算法給出了最優(yōu)答案。

尺度模型不僅解決了形態(tài)學(xué)中去噪的問題，而且解決了無監(jiān)督聚類問題。下面我們給出另外一種基于尺度模型的無監(jiān)督聚類算法，用于選擇合適的聚類數(shù)目。

2.3 基于互信息最大化的圖像分割

2.3.1 圖像分割算法的描述

絕大多數(shù)無監(jiān)督分割算法需要一些關(guān)于選擇NC方面的先驗(yàn)知識。如，我們事先知道腦圖像磁共振成像(MRI)可分為腦脊液、腦白質(zhì)、灰質(zhì)和背景。然而，沒有特定的專業(yè)領(lǐng)域知識，我們通常很難甚至不可能知道真正的聚類數(shù)目。

本文在這一節(jié)中給出了一種在無監(jiān)督圖像分割中解決選擇合適或者最優(yōu)NC問題的算法，其中水平集函數(shù)通常用來描述對比度不變的屬性[33]。一幅水平為t的圖像通常寫成，其中t∈R。這里，把一個(gè)特定的水平t當(dāng)作圖像分割的閥值。通過采用霍夫曼熵編碼方法，采用預(yù)先設(shè)定的準(zhǔn)則，算法可以把一幅圖像遞歸細(xì)分成相似的區(qū)域。該算法和它的分析是基于章節(jié)2.1描述的信息理論函數(shù)以及霍夫曼熵編碼方法。在每個(gè)階段對圖像直方圖進(jìn)行霍夫曼規(guī)則排序，保證了本文的互信息最大化的解不是局部的，而是全局最優(yōu)的。用i-Se標(biāo)記本文提出的基于互信息最大化的分割方法，它的偽代碼描述如下：

Algorithm:i-Se

Parameter:

Denote T as the number of maximum intensity bin for an image u(x), S as the clustering operation.

Do fort＝T—1,…,2,1

1. Calculate probability density

2. Sort pu(i) in descending order;

3. Combine the lowest two ordered probabilities, and produce a weighted intensity as follows:

Now the probability of intensity k becomes:

where

And the following diagram illustrates our clustering process based on Huffman coding strategy:

4. Calculate following mutual information:

2.3.2 關(guān)于圖像分割算法的分析

對于i-Se算法分析的主要目的是分析基于公式(25)的互信息函數(shù)單調(diào)性和光滑性，并進(jìn)一步給出對于 i-Se 算法的解是凸的。從以下兩個(gè)不等式開始(以下兩個(gè)引理的證明見附錄)。

引理1：令 p,q為概率密度函數(shù)，r＝p+q，那么有：—rlog(r)＜—(plog(p)+qlog(q))。

引理2：令 Stu(x)為圖像u(x)在水平t上的分割圖像，那么，聯(lián)合熵H(u,Stu)在任何聚類t水平都保存不變。

由這兩個(gè)引理，我們可以得出下面的命題。

命題3：令 Stu(x)為隨機(jī)變量u(x)在聚類水平t上的分割圖像，那么，隨著尺度變量t的衰減，I(u,Stu)的互信息是單調(diào)減少的。

證明：為了證明上面的命題，我們只需要證明互信息I(u,Stu)比H(u,St—1u)大。由(1)式有：I(u,Stu)＝H(u)+H(Stu)—H(u,Stu)，并且I(u,St—1u)＝H(u)+H(St—1u)—H(u,St—1u)。根據(jù)香濃熵的定義和引理1有：H(Stu)＞H(u,St—1u)。根據(jù)已被證明的引理2有：I(u,Stu)＞I(u,St—1u)。

在i-Se算法中，哈夫曼編碼樹的構(gòu)建對獲得全局最小解起著非常重要的作用。否則，公式(25)互信息中的凸性得不到保證。這種非凸的解通常存在于經(jīng)典的FCM算法，并且可以從電腦斷層掃描(CT)的體模圖像分類實(shí)驗(yàn)中觀察到。也即，F(xiàn)CM的解決方法是局部最小值。為了尋找全局的解決方法，i-Se算法采用類似于哈夫曼冗余編碼中的排序規(guī)則，通過合并兩個(gè)最小的概率密度值來構(gòu)建每一步最小的信息熵。同時(shí)，由命題3 可以得到：

公式(27)表明：(25)式的互信息計(jì)算存在一些冗余項(xiàng)，因此，可以通過減少冗余項(xiàng)的計(jì)算來增加i-Se算法的計(jì)算速度。

Algorithm:i-Se(Simpli fied version)

Parameters:

Denote T as the number of maximum intensity bin for an image u(x),

S as the clustering operation.

Calculate entropy Ht(Stu)at scale t＝T.

Do fort＝T—1,…,2,1

1. Calculate probability density

2. Sortpu(i) in descending order;

3. Combine the lowest two ordered probabilities and produce a weighted intensity as in (22);

4. Calculate Shannon entropy as follows:

3 實(shí)驗(yàn)結(jié)果

3.1 去噪結(jié)果

在多尺度i-De去噪算法的實(shí)驗(yàn)中主要關(guān)注具有非均勻光照背景的圖像，可以通過高帽變換有效地去除原始圖像存在的非均勻光照[46]。變換后的圖像被進(jìn)一步分割，以驗(yàn)證i-De尺度選擇算法。

圖1(a)是一幅帶有不均勻光照的原始米粒圖像。運(yùn)用閾值分割方法對原始圖像進(jìn)行分割，圖1(d)中給出了分割結(jié)果。從中可以看到，圖像底部的米粒不能理想地從非均勻背景中提取出來。Gonzalez等[46]曾指出只要結(jié)構(gòu)元素足夠大，對原始圖像進(jìn)行開運(yùn)算后，非均勻光照能夠被估算出來，但并沒有給出具體的尺度選擇方法。圖1(b)詮釋了尺度為6時(shí)的開操作圖像，在這個(gè)尺度公式(20)和(21)的互信息相交，可以看到，在這個(gè)交點(diǎn)的尺度不是最優(yōu)的，圖1(b)中用圓圈標(biāo)記出來的米粒的估計(jì)是不正確的。i-De尺度選擇算法選擇尺度為10，理想的結(jié)果見圖1(c)和圖1(f)。圖2表明隨著尺度t的增加，公式(20)和(21)的互信息將會(huì)相交。并且，它驗(yàn)證了在命題1中的I1和I2的預(yù)計(jì)：I1單調(diào)減少而 I2單調(diào)增加。然而，I1和I2不是凸函數(shù)，它的非凸性在圖3中被表示出來了，其中它詮釋了I1和 I2振蕩的差異。

圖1 米粒圖像去噪Fig.1. Denoising for rice grains image

圖2 各種不同尺度下米粒圖像的(20)和(21)的互信息Fig.2. Mutual information of (20) and (21) for grains image at various scales

圖3 圖2中互信息的差Fig.3. Difference of mutual information for Fig. 2

圖4為血管圖像去噪結(jié)果。其中，(a)～(d)中的Ttu對應(yīng)尺度分別為0、5、15、30；(e)～(h)為用水平集方法得到與(a)～(d)對應(yīng)的分割結(jié)果；(i)～(k)表明了非均勻背景的估算：Otu相應(yīng)的尺度t＝5，15，30。圖4(a)是具有非均勻背景的X射線血管醫(yī)學(xué)圖像。圖4(e)為用水平集方法[47]的分割結(jié)果。水平集在在血管左側(cè)底部停止演化，表明水平集方法不能有效地對具有非均勻背景的血管圖像進(jìn)行分割。本文用高帽變換對血管圖像中的非均勻背景進(jìn)行消除。公式(20)和(21)的互信息和它們之間的差異分別見圖5和圖6，i-De尺度分割算法預(yù)測最優(yōu)尺度為15，為了和其他尺度進(jìn)行對比，選擇具有代表性的尺度5和30。圖4中(a)～(d)分別為尺度0，5，15，30處的高帽變換圖像，(e)～(h)為用水平集法對濾波后的圖像分割結(jié)果，(i)～(k)顯示了對不均勻背景的估計(jì)。從這些圖像可以看出，如果選擇尺度過小，將會(huì)導(dǎo)致對非均勻背景的低估，正如圖4(f)圓圈標(biāo)記的那樣，它將會(huì)導(dǎo)致血管邊緣的泄露。另一方面，選擇尺度過大，將導(dǎo)致對背景的過高估計(jì)。只有選擇合適的尺度才會(huì)產(chǎn)生如4(g)那樣理想的分割結(jié)果。

圖4 血管圖像去噪Fig.4. Denoising for vessel image

3.2 分割結(jié)果

測試i-Se算法的第一個(gè)實(shí)驗(yàn)是在一個(gè)Shepp-Logan體模圖像上進(jìn)行的，它被實(shí)驗(yàn)者廣泛地用于驗(yàn)證它們?nèi)我舛SCT重建算法的準(zhǔn)確性。它具有六個(gè)已知的組織類型，其中的強(qiáng)度值分別為0，26，51，77，102和255。

圖7和圖8為用FCM和i-Se聚類方法得到的分割結(jié)果。在圖9和10為互信息和NC之間的關(guān)系，其中公式(25)和(26)之間的互信息被繪制了出來。它驗(yàn)證了定理3中的結(jié)論：I1的互信息是單調(diào)增長并隨著NC的增長匯聚于1。同時(shí)，I2單調(diào)減少并隨著 NC后達(dá)到0。

圖5 各種不同尺度下血管圖像的(20)和(21)互信息Fig.5. Mutual information of (20) and (21) for vessel image at various scales

圖6 圖5中互信息的差Fig.6. Difference of mutual information for Fig.5

圖7 基于FCM算法的體模圖像分割Fig.7. Segmentation of head phantom image using FCM method

圖8 基于i-S算法的體模圖像分割Fig.8. Segmentation of phantom image using i-Se algorithm

圖9 針對體模圖像FCM算法的互信息量Fig.9. Mutual information of FCM algorithm for head phantom image

圖10 針對體模圖像i-Se算法的互信息量Fig.10. Mutual information of i-Se algorithm for head phantom image

圖11 體模圖像中FCM和i-Se算法在互信息方面的比較Fig.11. Comparison of mutual information for FCM and i-Se algorithm

圖12 體模圖像中FCM和i-Se算法在互信息差方面的比較Fig.12. Comparison of the difference of mutual information for FCM andi-Sealgorithm

圖11表明當(dāng)NC是2或者3的時(shí)候i-Se算法的互信息量比FCM聚類算法的要多，這意味著i-Se算法的分割結(jié)果包含更多的圖像細(xì)節(jié)，這可以通過比較圖7(b)、(c)和圖8(b)、(c)得到。同時(shí)，F(xiàn)CM方法是一種局部算法，而且它會(huì)根據(jù)聚類中心不同的初始化產(chǎn)生不同的分割結(jié)果。這種局部現(xiàn)象在圖12中被很好地詮釋了，其中這種局部算法不能保證互信息差的單調(diào)性。相反，圖12表明了i-Se算法的互信息不同的曲線更平滑。更重要的是，它是單調(diào)減少的。

第二個(gè)測試i-Se算法的實(shí)驗(yàn)是在一幅真實(shí)的腦CT圖像上進(jìn)行的。它不是一幅合成圖像并且我們沒有預(yù)先知道準(zhǔn)確的NC。圖13為測試的CT圖像及用FCM和i-Se算法得到的分割結(jié)果?；バ畔⒑蚇C之間的關(guān)系以及不同分別在圖14和圖15中表示了出來。如圖14所示，隨著NC的增加，i-Se和FCM算法的等式(25)的互信息都在增加。但是也可以觀察出來，在NC相同的情況下，i-Se算法的互信息量比FCM算法的大。同時(shí)，與FCM的分割結(jié)果圖像相比，更詳細(xì)的細(xì)節(jié)信息被保存在i-Se分割圖像中。隨著NC的增加，這種優(yōu)勢更加明顯。圖15中，隨著NC 的增加，F(xiàn)CM算法的互信息的差異的是減少振蕩的，而隨著NC的減少i-Se算法單調(diào)減少。

i-Se算法的第三個(gè)實(shí)驗(yàn)是在Lena圖上進(jìn)行的。實(shí)驗(yàn)結(jié)果見圖16～18。把腦T圖像和Lena圖像之間的結(jié)果進(jìn)行比較，我們的方法可以有效地保存分割圖像的詳細(xì)信息。在其它測試圖像上進(jìn)行了更多的實(shí)驗(yàn)，用來驗(yàn)證i-Se算法在揭示原始圖像和分割圖像內(nèi)在關(guān)系上的有效性。

圖13 CT圖像分割Fig. 13. Segmentation for CT image

圖14 CT圖像中 FCM 和 i-Se 算法在互信息方面的比較Fig.14. Comparison of mutual information for CT image using FCM and i-Se algorithm

圖15 CT圖像中FCM和i-Se算法在互信息差方面的比較Fig.15. Comparison of the difference of mutual information for CT image using FCM and i-Sealgorithm

圖16 Lena圖像分割Fig.16. Segmentation for Lena image

圖17 Lena 圖像中 FCM 和 i-Se 算法在互信息方面的比較Fig.17. Comparison of mutual information for Lena image using FCM and i-Se algorithm.

圖18 Lena 圖像中 FCM 和 i-Se 算法在互信息差方面的比較Fig.18. Comparison of the difference of mutual information for Lena image using FCM and i-Se algorithm

4 結(jié)果和討論

在本文中，我們利用互信息作為在多尺度分析中最優(yōu)尺度選擇的相似性的衡量標(biāo)準(zhǔn)，應(yīng)用到圖像去噪和分割中。首先，我們表示了在多尺度分析中的圖像處理算子并提出了互信息函數(shù)。對于形態(tài)學(xué)開操作算子，我們給出了一個(gè)i-De尺度選擇算法，它的好處是使在多尺度表示的圖像中讓多尺度的選擇自動(dòng)化。對于聚類算子，我們給出了一種i-Se聚類算法，該算法通過把一個(gè)特定水平t當(dāng)做圖像分割的閥值。對于i-De和i-Se算法的可行性，我們對目標(biāo)函數(shù)的單調(diào)性和平滑性做了詳細(xì)的數(shù)學(xué)分析。

利用i-De算法，通過選擇合適的尺度來消除不均勻背景。這種非均勻的光照在大多數(shù)相機(jī)或醫(yī)療圖像中是一種常見的現(xiàn)象。i-Se算法輸出了聚類分割的合適水平，在這樣的尺度下，我們可以得到分割圖像的最大量信息。

i-Se算法的一個(gè)缺點(diǎn)是模型中只考慮到灰色特征，其他的特征(如：像素之間的空間相關(guān)性、邊界和紋理)沒有被考慮到，需要一種更詳細(xì)的局部模型來產(chǎn)生更可靠的結(jié)果。

許多其他的建模任務(wù)，如：化學(xué)、物理和生物建模，可以在尺度空間建模方面進(jìn)行探索。然而，本研究的主要目標(biāo)是為圖像去噪和分割任務(wù)建立一種尺度選擇標(biāo)準(zhǔn)。在未來，我們將會(huì)就關(guān)于我們的算法應(yīng)用到其他的模型任務(wù)中做更多的工作，并且針對模型中現(xiàn)存的大多數(shù)限制進(jìn)行研究解決。

附錄

引理1：令p,q為概率密度函數(shù)，r＝p+q，那么有：—rlog(r)＜—(plog(p)+qlog(q))。

證明：因?yàn)閜,q為概率函數(shù)，所以有和對于一個(gè)遞增函數(shù)f(x)＝xa,a∈[0,1],x∈(0,∞)，可以得到：(p+q)p＞pp和(p+q)q＞pq。所以有：rr＝(p+q)p+q＝(p+q)p*(p+q)q＞pp*qq。兩邊取對數(shù)可得：rlog(r)＞plog(p)+qlog(q)。

引理2：令 Stu(x)為圖像u(x)在水平t上的分割圖像，那么，聯(lián)合熵H(u,Stu)在任何聚類t水平都保存不變。

證明：對于自聯(lián)合熵H(u,u)，我們有puu(i,j)i＝j(luò)≠0和puu(j,j)i≠j＝0。即自聯(lián)合熵是一個(gè)對角矩陣。假設(shè)合并u(x)任意兩個(gè)概率密度，如i和j，那么會(huì)產(chǎn)生一個(gè)如公式(22)的加權(quán)密度函數(shù)，以及Stu(x)的概率密度k變成了公式(23)。所以在對角矩陣H(u,u)上的聯(lián)合直方圖有如下的對角元素移動(dòng)：puu(k,k)→puStu(k,k)，puu(i,i)→puStu(i,k)和puu(j,j)→puStu(j,k)。由于聚類的結(jié)果只在聯(lián)合分布上移動(dòng)元素，并不改變聯(lián)合熵的整個(gè)值。即有H(u,u)＝H(u,Stu)。相似地，有H(u,u)＝H(u,St—1u)。因此有H(u,Stu)＝H(u,St—1u)，從而該引理得證。

[1] Spohn H. Large Scale Dynamics of Interacting Particles [M]. New York: Springer Berlin Heidelberg, 1991.

[2] Katsoulakis MA, Majda AJ, Vlachos DG. Coarsegrained stochastic processes and Monte Carlosimulations in lattice systems [J]. Journal of Computational Physics, 2003, 186: 250-278.

[3] Xiao L, Chen LX, Xiao JZ. A new algorithm for shortest path problem in large-scale graph [J]. Applied Mathematics Information Sciences, 2012, 6(3): 397-400.

[4] Pappalardo F, Palladini A, Pennisi M, et al. Mathematical and computational models in tumor immunology [J]. Mathematical Modelling of Natural Phenomena, 2012, 7(3): 186-203.

[5] Keunings R. Micro-macro methods for the multiscale simulation of viscoelastic flows using molecular methods of kinetic theory [J]. Rheology Reviews, 2004, 1: 67-98.

[6] Erban R, Othmer HG. From individual to collective behavior in bacterial chemotaxis [J]. SIAM Journal on Applied Mathematics, 2004, 65(2): 361-391.

[7] Alemani D, Pappalardo F, Pennisi M, et al. Combining cellular automata and lattice Boltzmann method to model multiscale avascular tumor growth coupled with nutrient diffusion and immune competition [J]. Journal of Immunological Methods, 2012, 376(1-2): 55-68.

[8] Bradley C, Bowery A, Britten R, et al. Open CMISS: a multi-physics & multi-scale computational infrastructure for the VPH/Physiome project [J]. Progress in Biophysics and Molecular Biology, 2011, 101(1): 32-47.

[9] Bellouquid A, Bianca C. Modelling aggregationfragmentation phenomena from kinetic to macroscopic scales [J]. Mathematical and Computer Modelling, 2010, 52(5-6): 802-813.

[10] Marr D. Vision [M]. Dallas: Freeman Publishers, 1982.

[11] Koenderink JJ. The structure of images [J]. Biological Cybernetics, 1984, 50(5): 363-370.

[12] Lindeberg T. Scale-space Theory in Computer Vision [M]. Norwell: Kluwer Academic Publishers, 1994.

[13] Witkin AP. Scale-space fi ltering [C] // Proceedings of the 8th International Joint Conference on Artif i cial Intelligence, 1983, 2: 1019-1022.

[14] Babaud J, Witkin AP, Baudin M, et al. Uniqueness of the Gaussian kernel for scale-space fi ltering [J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 1986, 8: 26-33.

[15] Mallat S. A Wavelet Tour of Signal Processing [M]. Utah: Academic Press, 1988.

[16] Portilla J, Strela V, Wainwright MJ, et al. Image denoising using scale mixtures of Gaussians in the wavelet domain [J]. IEEE Transactions on Image Processing, 2003, 12(11): 1338-1351.

[17] Lewis AS, Knowles G. Image compression using the 2-D wavelet transform [J]. IEEE Transactions on Image Processing, 1992, 1(2): 244-250.

[18] Adami N, Signoroni A, Leonardi R. State-of-theart and trends in scalable video compression with wavelet-based approaches [J]. IEEE Transactions on Circuits and Systems for Video Technology, 2007, 17(9): 1238-1255.

[19] Allen RL, Kamangar FA, Stokely EM. Laplacian and orthogonal wavelet pyramid decompositions in coarse-to-f i ne registration [J]. IEEE Transactions on Signal Processing, 1993, 41(12): 3536-3541.

[20] Sapiro G. Geometric Partial Differential Equations and Image Analysis [M]. Cambridge: Cambridge University Press, 2001.

[21] Perona P, Malik J. Scale-space and edge-detection using anisotropic diffusion [J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 1990, 12(7): 629-639.

[22] Rudin LI, Osher S, Fatemi E. Nonlinear total variation based noise removal algorithms [C] // Proceedings of the 11th Annual International Conference of the Center for Nonlinear Studies on Experimental mathematics : Computational Issues in Nonlinear Science, 1992: 259-268.

[23] Black MJ, Sapiro G, Marimont DH, et al. Robust anisotropic diffusion [J]. IEEE Transactions on Image Processing, 1998, 7(3): 421-432.

[24] Tschumperle D. Fast anisotropic smoothing of multi-valued images using curvature-preserving PDE’s [J]. International Journal of Computer Vision, 2006, 68(1): 65-82.

[25] Vincken KL, Koster ASE, Viergever MA. Probabilistic multiscale image segmentation [J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 1997, 19(2): 109-120.

[26] Leung Y, Zhang JS, Xu ZB. Clustering by scalespace filtering [J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2000, 22(12): 1396-1410.

[27] Petrovic A, Escoda OD, Vandergheynst P. Multiresolution segmentation of natural images: from linear to nonlinear scale-space representations [J]. IEEE Transactions on Image Processing, 2004, 13(8): 1104-1114.

[28] Bertalmio M, Sapiro G, Caselles V, et al. Image inpainting [C] // Proceedings of SIGGRAPH 2000, 2000: 417-424.

[29] Bertalmio M, Bertozzi AL, Sapiro G. Navier-stokes, fl uid dynamics, and image and video inpainting [C] // Proceedings of the 2001 IEEE Computer Society Conference on Computer Vision and Pattern Recognition, 2001, 1: 355-362.

[30] Chan TF, Shen JH. Variational image inpainting [J]. Communications on Pure and Applied Mathematics, 2005, 58(5): 579-619.

[31] Serra J. Image Analysis and Mathematical Morphology [M]. Utah: Academic Press, 1982.

[32] Jackway PT, Deriche M. Scale-space properties of the multiscale morphological dilation erosion [J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 1996, 18(1): 38-51.

[33] Guichard F, Moisan L, Morel JM. A review of P.D.E. models in image processing and image analysis [J]. Journal de Physique IV France, 2002, 12: 137-154.

[34] Lindeberg T. Edge detection and ridge detection with automatic scale selection [C] // Proceedings of IEEE Computer Society Conference on Computer Vision and Pattern Recognition, 1996: 465-470.

[35] Dolcetta IC, Ferretti R. Optimal stopping time formulation of adaptive image filtering [J]. Applied Mathematics and Optimization, 2001, 43(3): 245-258.

[36] Solo V. A fast automatic stopping criterion for anisotropic diffusion [C] // IEEE International Conference on Acoustics, Speech, and Signal Processing, 2002, 2: 1661-1664.

[37] Mrazek P, Navara M. Selection of optimal stopping time for nonlinear diffusion filtering [J]. International Journal of Computer Vision, 2003, 52(2-3): 189-203.

[38] Papandreou G, Maragos P. Image denoising in nonlinear scale-spaces: automatic scale selection via cross-validation [C] // IEEE International Conference on Image Processing, 2005: 1033-1036.

[39] Papandreou G, Maragos P. A cross-validatory statistical approach to scale selection for image denoising by nonlinear diffusion [C] // IEEE Computer Society Conference on Computer Vision and Pattern Recognition, 2005: 625-630.

[40] Sun J, Xu ZB. Scale selection for anisotropic diffusion fi lter by Markov random fi eld model [J]. Pattern Recognition, 2010, 43(8): 2630-2645.

[41] Wen TX, Gu J, Zhang ZQ, et al. Scale selection for morphological top-hat transformation based on mutual information [C] // The 3rd International Congress on Image and Signal Processing, 2010(3): 1092-1096.

[42] Rigau J, Feixas M, Sbert M, et al. Medical image segmentation based on mutual information maximization [C] // MICCAI Lecture Notes in Computer Science, 2004, 3216: 135-142.

[43] Lu QW, Chen WF. Unsupervised segmentation of medical image based on difference of mutual information [J]. Science in China, Series F: Information Sciences, 2006, 49(4): 484-493.

[44] Lu ZT, Feng QJ, Shi PC, et al. Unsupervised segmentation of medical image based on FCM and mutual information [C] // IEEE International Conference on Complex Medical Engineering, 2007, 1: 513-516.

[45] Bezdek JC. Pattern Recognition with Fuzzy Objective Function Algorithms [M]. New York: Plenum Press, 1981.

[46] Gonzalez RC, Woods RE, Eddins SL. Digital Image Processing Using MATLAB [M]. Beijing: Publishing House of Electronics Industry, 2004.

[47] Chan T, Vese L. An active contour model without edges [J]. IEEE Transactions on Image Processing, 2001, 10: 266-277.

Image Denoising and Segmentation Based on Mutual Information Criterion

WEN Tiexiang1,2PAN Zhengyang1,3GU Jia1
1( Shenzhen Institutes of Advanced Technology, Chinese Academy of Sciences, Shenzhen 518055, China )

2( University of Chinese Academy of Science, Beijing 100049, China )

3( University of Science and Technology of China, Hefei 230000, China )

Scalespace play an important role in many computer vision tasks. Automatic scale selection is the foundation of multi-scale image analysis, but its performance is still very subjective and empirical. To automatically select the appropriate scale for a particular application, a scale selection model based on information theory was proposed in this paper. The proposed model utilizes the mutual information as a measuring criterion of similarity for the optimal scale selection in multi-scale analysis, with applications to the image denoising and segmentation. Firstly, the multi-scale image smoothing and denoising method based on the morphological operator was studied. This technique does not require the prior knowledge of the noise variance and can effectively eliminate the changes of illumination. Secondly, a clusteringbased unsupervised image segmentation algorithm was developed by recursively pruning the Huffman coding tree. The proposed clustering algorithm can preserve the maximum amount of information at a speci fi c clustering number from the information-theoretical point of view. Finally, for the feasibility of the proposed algorithms, its theoretical properties were analyzed mathematically and its performance was tested through a series of experiments, which demonstrate that it yields the optimal scale for the developed image denoising and segmentation algorithms.

scale selection; mutual information; denoising; segmentation

TP 391

A

2013-11-20

溫鐵祥，博士研究生，研究方向?yàn)獒t(yī)學(xué)圖像處理和計(jì)算機(jī)視覺；潘正洋，碩士研究生，研究方向?yàn)獒t(yī)學(xué)圖像處理；辜嘉(通訊作者)，研究員，博士生導(dǎo)師，研究方向?yàn)橛?jì)算機(jī)視覺和醫(yī)學(xué)成像技術(shù)，E-mail：jia.gu@siat.ac.cn。