宋金濤，于 晗，潘振寬

(1. 青島大學(xué)計算機科學(xué)技術(shù)學(xué)院，山東 青島 266071；2. 國防科技大學(xué)計算機學(xué)院，湖南 長沙 410000)

0 引言

通過對設(shè)計的能量泛函取極值, 可有效實現(xiàn)圖像恢復(fù)[1-2]，包括圖像噪聲去除、 圖像修復(fù)、 圖像去模糊等. 相關(guān)能量泛函通常包含數(shù)據(jù)項和規(guī)則項，前者表示恢復(fù)后的圖像與原圖像的逼近程度，后者是恢復(fù)后圖像特征保持的關(guān)鍵. 經(jīng)典的圖像恢復(fù)變分模型的規(guī)則項基于一階導(dǎo)數(shù)，但在結(jié)果中會呈現(xiàn)階梯效應(yīng). 為克服一階導(dǎo)數(shù)模型的缺陷，基于二階導(dǎo)數(shù)規(guī)則項變分圖像恢復(fù)模型的研究成為近年國際圖像處理與分析領(lǐng)域研究的熱點之一. 但由于該類模型會出現(xiàn)四階偏微分方程，其計算效率成為該類方法實際應(yīng)用的主要瓶頸[3].

定義于圖像空間中的函數(shù)，其二階導(dǎo)數(shù)具有多種形式，它們在變分圖像恢復(fù)模型中表現(xiàn)出的特點亦有差異. 該領(lǐng)域已涌現(xiàn)出基于拉普拉斯[4]、 曲率[5]、 海森矩陣[6]、 對稱海森矩陣[7]、 歐拉彈性能(Euler’s elastica)的模型[8]，并應(yīng)用于圖像噪聲去除、 圖像修復(fù)、 圖像分割、 虛幻輪廓恢復(fù)等問題[3]. 本研究擬針對變分圖像恢復(fù)模型僅引入1個對偶變量和1個約束方程，并結(jié)合簡單的約束投影方法，提出一種快速對偶投影算法.

1 相關(guān)研究

包含歐拉彈性能的能量項, 最早由Nitzberg等[9]應(yīng)用于深度圖像分割與虛幻輪廓恢復(fù)；Masnou等[10]將其推廣到大破損圖像修復(fù)；Zhu等[11]將其推廣到圖像分割變分模型，其形式為

(1)

式中：Ω為有界圖像區(qū)域；a、b為兩項平衡的懲罰參數(shù)；u為恢復(fù)后圖像強度；對于平面圖像，dx=dx1dx2.

歐拉彈性能的非光滑性、 非線性、 非凸性及其高階導(dǎo)數(shù)等特點導(dǎo)致了上述模型計算上的困難，主要體現(xiàn)在由于變分導(dǎo)致的高階非線性偏微分方程的離散復(fù)雜、 迭代計算的低效和不穩(wěn)定性等方面. 上述不同模型研究的同時，先后出現(xiàn)了梯度降方法[8, 12]、 圖割方法[13]、 增廣拉格朗日方法(ALM)[14]、 對偶ALM方法[15]、 算子分裂方法[16]等. Kang等[3]以圖像分割模型為背景對上述算法進行了比較系統(tǒng)的回顧. 早期的梯度降方法[8]需要求解強非線性耦合的四階偏微分方程，盡管Chan等[8]采用預(yù)置條件方法，文獻[12]采用離散梯度方法改善計算的穩(wěn)定性，但計算效率仍較低. Tai等[14]通過引入多個輔助變量和拉格朗日乘子將原模型轉(zhuǎn)化為一系列可簡單求解的子問題交替求解，其計算框架為增廣拉格朗日方法(ALM)，或稱交替方向乘子方法(ADMM). 變量的增加不僅加重了存儲負擔(dān)，還增加了許多額外的計算量. 為減少引入的輔助變量數(shù)目，Zhang等[15]提出了基于對偶方法和ALM方法的混合方法，與Tai等[14]的方法相比僅減少了1個變量和1個乘子. 最新的算子分裂方法[16]通過引入多個指示函數(shù)對原模型進行修改，并對新的模型采用算子分裂方法多步求解，具有較高效率，但其松弛技術(shù)將對計算結(jié)果產(chǎn)生影響. 圖割方法[13]采用各向異性模型對原模型近似，盡管計算效率高，但得到的并非原模型的解.

對于給定含噪聲圖像f(x)∈Ω，x∈Ω，基于歐拉彈性能的圖像恢復(fù)變分模型[8]為

(2)

其中，右端第二項為數(shù)據(jù)項，用于刻劃恢復(fù)后清晰圖像u與原圖像f的接近程度.

為了克服光滑區(qū)域梯度過小引起的計算穩(wěn)定性問題，Chan等[8]對變分法求得的圖像擴散的梯度降方程進行預(yù)置條件處理，將其修改為

(3)

從而可采用迎風(fēng)差分格式對其離散求解，并可選取較大的時間步長，但計算效率依然較低.

為了克服上述問題，Tai等[14]通過引入多個輔助變量、 拉格朗日乘子和懲罰參數(shù)，將式(2)轉(zhuǎn)化為交替優(yōu)化方法，即增廣拉格朗日方法(ALM).

引入輔助變量w=·n, |m|≤1，則式(2)轉(zhuǎn)化為受約束的多變量優(yōu)化問題

(4)

s.t.w-·n=0，m-n=0，|m|≤1

(5)

但純ALM方法引入過多的輔助變量和拉格朗日乘子. 為了降低求解規(guī)模，Zhang等[15]首先采用Legendre-Fenchel變換將歐拉彈性能轉(zhuǎn)化為

(6)

并進一步引入約束C(u，n)=(n?n-I), 將式(2)轉(zhuǎn)化為更加簡單的能量泛函極值的迭代優(yōu)化形式. 從而比Tai等[14]的方法減少了1個優(yōu)化變量和3個拉格朗日乘子. 盡管該方法相關(guān)子優(yōu)化問題比Tai等[14]提出的方法復(fù)雜，總體計算效率較高.

本研究設(shè)計的對偶投影方法是圖像恢復(fù)經(jīng)典總變差(total variation，TV)模型方法[17]的推廣，該方法由Chambolle[18]提出. 其基本思路是將TV模型

(7)

轉(zhuǎn)化為如下鞍點問題

(8)

其中:p為對偶變量，源于總變差的對偶表達形式

(9)

固定p，求關(guān)于u的變分問題，得

u-f-γ·p=0

(10)

將式(10)帶入式(8)，求關(guān)于p的極大值問題

(11)

依據(jù)受不等式約束優(yōu)化問題的KKT條件、 梯度降方法及半隱式迭代可得到關(guān)于pk+1的迭代計算格式.

2 含歐拉彈性能的圖像恢復(fù)模型的對偶投影方法

為了應(yīng)用對偶投影法求解復(fù)雜的變分模型(2)，首先將其轉(zhuǎn)化為如下迭代優(yōu)化形式

(12)

其中

(13)

采用對偶方法，式(12)可轉(zhuǎn)化為

(14)

固定p，求關(guān)于u的極小值問題，得到公式(10). 固定u，求關(guān)于p的極大值問題，其相關(guān)梯度升方程

(15)

采用一階差分離散(15)可得

(16)

其中:pk+1，0=pk；gk=1，0=gk.

式(16)的投影公式為

(17)

直到

(18)

ε為誤差容限，本研究取ε=0.001，τ≤0.125.

令qi，j=(則wi，j=·ni，j，由此可得

基于上述離散，迭代公式(17)可簡潔表達為

(19)

3 數(shù)值實驗

實驗的目的是對提出的對偶投影方法(dual projection, DP)的計算效率進行比較驗證，擬比較的方法為經(jīng)典的梯度降方法(GD)[8]，增廣拉格朗日方法(THC)[14]，及快速對偶增廣拉格朗日方法(augmented lagrangian primal dual，ALPD)[15]. 但基于Tai等[14]的研究，THC方法的效率是GD算法的近500倍，而Zhang等[15]表明ALPD算法計算效率大約是THC的3倍.

在簡單圖像和真實圖像噪聲去除問題中應(yīng)用本算法，詳見圖1、 2所示. 所有的圖像都是灰度的，運行環(huán)境為Matlab 2016a, 參數(shù)設(shè)置均為a=0.5，b=1，γ=1. 在圖1、 2中，原始圖像受到高斯白噪聲的影響，均方差為0.001. 將DP算法與THC算法進行對比發(fā)現(xiàn)，除了在處理復(fù)雜圖像時的邊界略有差異外，兩種算法都具有良好的去噪效果，而DP算法在邊界處理和對比度保持方面有更好的細節(jié).

圖1 兩幅簡單圖像的去噪結(jié)果 Fig.1 Denoising results of two simple images

圖2 兩幅真實圖像的去噪結(jié)果 Fig.2 Denoising results for two real images

表1為兩種算法對圖1中2幅圖像去噪的效果 (通過信噪比PSNR)和效率(通過迭代次數(shù)以及運行所需的秒數(shù))的比較. 每幅圖像都經(jīng)過4種尺寸的測試，以更準確地比較差異. 結(jié)果表明，該方法需要較少的時間來實現(xiàn)相同甚至更好的去噪質(zhì)量，而且在實驗對象為256 px×256 px的“Lena”圖像時，每一步迭代 THC算法需要0.074 s，但DP算法只需要0.012 s. 主要原因是該算法只有u和p需要迭代求解. 從整體上看，DP算法在效率和結(jié)果上都有較大的提高.

針對圖2中2幅真實圖像實驗結(jié)果以相同的方式列入表2，獲得的結(jié)論與對圖1的結(jié)論一致.

表1 DP算法和THC算法對不同簡單圖像處理性能比較

表2 DP算法和THC算法對不同真實圖像處理性能比較

圖3給出了大小為256 px×256 px的“Lena”圖像去噪的能量變化和信噪比. 其中圖3(a)為THC算法對 “Lena”圖像的能量變化與迭代步數(shù)的關(guān)系曲線；圖3(b)為DP算法的相應(yīng)曲線；兩種算法恢復(fù)出清晰圖像所需迭代次數(shù)差異微小，但每步計算時間差別較大. 圖3(c)為THC算法與DP算法的峰值信噪比(PSNR)隨迭代次數(shù)變化曲線，可以看出，DP算法取得更高的峰值信噪比PSNR值.

圖3 THC算法與DP算法的能量函數(shù)以及信噪比圖Fig.3 THC algorithm and our DP algorithm energy function and signal to noise ratio diagram

圖4 不同噪聲環(huán)境下的去噪結(jié)果 Fig.4 Denoiseing results for different noise

圖4展示了“Lena”圖像在不同噪聲環(huán)境下的去噪結(jié)果，其中圖4(b)為噪聲密度為0.01的椒鹽噪聲，圖4(d)、 (e)、 (f)分別為DP算法的泊松噪聲、 椒鹽噪聲和斑點噪聲的去噪結(jié)果. 由結(jié)果可知，DP算法在不同的噪聲環(huán)境下均有良好的去噪結(jié)果，但是因為離散時采取的網(wǎng)格間距為1，所以不能很好地處理大塊的椒鹽噪聲.

實驗表明，在保持較高信噪比的情況下，DP方法計算效率大約是THC方法的5倍，即是ALPD算法的1.5倍. 盡管上述實驗數(shù)量有限，但在針對其他20余篇圖像噪聲去除的實驗中取得了相同的結(jié)論.

4 結(jié)語

研究針對含歐拉彈性能的噪聲去除變分圖像恢復(fù)模型, 提出一種快速的對偶投影算法，該方法僅計算2個變量，包含1個約束，且約束形式簡單，可用于投影方法的快速實現(xiàn)，通過數(shù)值實驗與現(xiàn)有快速算法比較，大大提高了計算效率. 該算法可直接推廣到含歐拉彈性能的圖像分割、 圖像深度分割、 虛幻輪廓恢復(fù)、 圖像修復(fù)等諸多問題求解.