1. 引言
随着医学成像技术的发展,医学图像在临床诊断和治疗规划中的应用越来越广泛。由于设备精度和环境因素等影响而引起的噪声和模糊,往往导致图像造成不同程度的退化。噪声和模糊的存在会导致图像细节丢失和质量下降,对医学图像的分析和诊断造成不利影响。因此,如何有效地从图像中去除模糊和噪声已经成为图像处理领域的重要研究课题,对医学图像质量的提升和推动医学影像技术的进步有着重要意义。模糊过程可以用点扩散函数(PSF)来表征,如果PSF已知,那么从带有模糊和噪声的图像中复原图像称为非盲去模糊,如果PSF未知则是盲去模糊。
去模糊是高度不适定的逆问题,可以通过应用先验信息并构建有效的正则项来解决图像复原难题。全变分算法(TV)是应用最为广泛的图像复原技术之一,Rudin等人[1]通过最小化正则化项,重建图像的边缘信息,从而去除模糊和噪声并提升了复原图像质量,很多图像复原方法都是在此种算法基础上发展而来。Richardson [2]和Lucy [3]考虑图像模糊的情况下对泊松噪声图像进行恢复,提出了一种基于对数似然函数最小化的去噪算法(RLTV),该方法因引入了泊松统计量,所以在泊松图像复原中得到了广泛的应用。由于反卷积问题的病态性,该算法在经过几次迭代后会放大噪声,造成阶梯伪影的情况。文献[4]将全变分重叠组稀疏作为已知先验来解决泊松噪声图像的去模糊问题,该方法在去除噪声和保留边缘的同时,减少了全变分正则化处理图像中普遍存在的阶梯效应。
利用小波变换后高频图像具有稀疏性的特征,许多研究人员将小波变换与全变分相结合进行图像复原,例如在文献[5]中,作者将双树复小波变换和
范数作为稀疏约束项,提出了两种修正正则化参数的估计方法,使参数的设置变得简单。类似的,Zhang等人[6]将小波框架和
范数作为正则化项进行图像复原,但是
范数是非凸问题,求解起来比较麻烦。F. Luisier等人[7] [8] [9]提出了非迭代的去模糊方法用于复原含有泊松噪声的模糊图像,该方法将反卷积过程线性化为初等函数的线性组合,每个初等函数由维纳滤波和小波变换阈值组成,最后用无偏估计来优化线性组合的系数。该方法的优点是在求解过程中不需要迭代,计算复杂程度较低。
总广义变分(TGV) [10] [11]和分数阶全变分(FOTV) [12]是全变分的改进算法,两者的性能相当,但由于TGV比FOTV涉及更多的未知变量,所以FOTV比TGV效率更高。在基于FOTV算法基础之上,Yifei [13]等人为了去除泊松噪声并保持图像的高阶平滑性,提出了基于分数阶全变分正则化和泊松噪声统计特性的反卷积模型,应用交替方向乘子法(ADMM) [14]和最大期望算法(EM) [15]求解模型得到了具有理论保证的高效算法。该方法虽然能在一定程度上保持图像的高阶平滑性,但去除不同程度的泊松噪声时,复原效果并不理想。针对这种情况,本文结合小波变换以保持底层图像的高阶平滑,提出了基于分数阶全变分正则化的去模糊方法,并建立了非盲去模糊模型,针对不同程度的泊松噪声得到了很好的复原效果。
2. 相关工作
图像模糊过程可以用线性方程来建模:
(1)
其中,u是清晰图像,H是模糊核,g是模糊图像,
表示二维卷积算子,本文针对非盲去模糊问题进行研究。
2.1. 泊松噪声
泊松噪声可以用泊松分布来描述,假设一个随机变量f遵循参数为
的泊松分布,概率密度函数由下式定义:
(2)
其中,参数
决定了f的期望值以及它的方差。假设在每个像素处的测量数据
遵循泊松分布,其基础真值是
,因此概率密度函数为:
(3)
其中,N是像素的总数,
,每个
,
,由贝叶斯定理可知,对于给定的f和g的后验概率密度为:
(4)
然后取(4)的负对数似然,我们得到
(5)
我们可以通过最大后验概率(MAP)来估计清晰图像u,即
。由于模糊可以表示为
,H表示模糊核是一个实矩阵,则后验概率
和先验概率
分别与
和
相关。利用对数似然形式,MAP估计可表示为:
(6)
其中,
表示图像域
离散网格上的指示函数,
定义为内积,
与
相关,称为图像先验。
2.2. FOTV正则化模型
本文将分数阶全变分视为图像先验,给定一个图像域
,将其离散为一个矩形网格
,那么图像可以表示为欧几里得空间
中的矩阵,记为
。给定分数阶数导数
,Zhang,J提出了分数阶导数[16],将离散分数阶梯度定义为:
(7)
其中沿x轴和y轴的分数阶导数
近似为:
(8)
(9)
其中,K是在每个像素处近似分数阶导数的相邻像素的个数。当
为Gamma函数,系数
定义为
,将u的离散分数阶全变分定义为:
(10)
其中,
表示1范数,根据
,对于
,Zhang,J [16]给出离散分数阶散度系数
:
(11)
将FOTV正则化作为图像先验
合并到MAP估计(6)中,当数据被泊松噪声破坏时,文献[12]提出以下FOTV正则化模型:
(12)
其中,h表示模糊核,
是一个加权参数,在正则化项和数据拟合之间取得平衡。
3. 本文模型及求解
本文提出融合小波变换的分数阶全变分正则化模型:
(13)
其中,h表示模糊核,W表示小波变换,
表示平衡参数。
本文采用交替方向乘子法(ADMM) [14]来求解所提出的模型(13)。引入以下三个辅助变量
,
和
,并将(13)式等价的表示为如下的(14)式:
(14)
相应的增广拉格朗日泛函由下列式子表达:
(15)
其中,
,
和
是对偶变量或拉格朗日乘子,
和
是三个正参数,
表示Frobenius范数。然后ADMM产生以下迭代:
(16)
为了解决u子问题,我们计算(15)相对于u的梯度,得到最优性条件:
(17)
这里
是
的伴随,它是(11)中定义的散度算子,
是h顺时针旋转90◦得到的h的伴随核,
是W的转置。在周期边界条件下,由于复合算子
、
和
各自的循环变换矩阵,可以用快速傅里叶变换(FFT)对角化。所以,我们给出的u的封闭形式解如下:
(18)
其中,,F表示傅里叶变换,
表示傅里叶逆变换,矩阵平方、除法和绝对值都是按分量执行的。收缩算子可给出的z子问题的封闭形式解:
(19)
其中,
,
表示分量乘法,所有算术运算符都以分量方式执行。求函数L对g的导数,我们得到了g子问题的最优性条件:
(20)
由于
,我们将上面的方程改写为二次方程:
(21)
并选择正解,即:
(22)
对于求解e子问题,可以得到下式:
(23)
对上式进行简化得到:
(24)
针对求解e子问题我们选取硬阈值迭代算法进行处理,下面先给出硬阈值迭代算子的表达式:
(25)
根据硬阈值算子,给出e子问题的封闭形式解如下:
(26)
本文算法流程在算法1中表述。
算法1. 融合小波变换的分数阶全变分去模糊算法
输入:
|
设置参数:
|
初始化:
|
当
时 |
通过(18)式子求解
通过(19)式子求解
通过(22)式子求解
通过(26)式子求解
输出:
|
4. 数值实验
本文选取MRI数据集进行了数值实验来证明所提出的方法在去除泊松噪声和高斯模糊伪影方面的有效性。本文用峰值信噪比(PSNR)和结构相似指数(SSIM)来评价图像去模糊的质量,PSNR定义为:
(27)
SSIM定义为:
(28)
其中,u代表原始图像,g代表去模糊后的复原图像,MES是u与g之间的均方误差,
是协方差,
是u的平均值,
是g的平均值,
是u的方差,
是g的协方差。
,
是维持稳定的常数,
,
,L是像素的动态范围。可以看出PSNR值越大,图像复原质量越好,SSIM值越大说明复原图像与原图差异性越小,复原效果越接近原图。
接下来对医学图像进行去模糊处理,并与RLTV [2] [3]、FOTV [12]、BM3D [17]三种去模糊方法进行了比较。
本文随机选取了MRI数据集的六张512 × 512大小人脑切片图像(Image A, Image B, Image C, Image D, Image E, Image F)进行数值实验,利用核大小为3 × 3和标准差为1的高斯模糊核以及峰值水平(peak)分别为200、255和300的泊松噪声对清晰图像进行退化。每个峰值对应一个级别的泊松噪声,峰值越小,图像看起来就越嘈杂,去模糊过程越具有挑战性。
Table 1. PSNR values recovered by different algorithms under different levels of Poisson noise
表1. 不同程度泊松噪声下不同算法复原的PSNR值
|
Peak = 200 |
Peak = 255 |
Peak = 300 |
Image |
Ours |
RLTV |
FOTV |
BM3D |
Ours |
RLTV |
FOTV |
BM3D |
Ours |
RLTV |
FOTV |
BM3D |
A |
35.51 |
34.64 |
35.44 |
34.98 |
35.91 |
34.98 |
35.80 |
35.24 |
36.21 |
35.18 |
36.17 |
35.42 |
B |
35.51 |
34.40 |
35.43 |
34.83 |
35.96 |
34.71 |
35.89 |
35.05 |
36.27 |
34.9 |
36.24 |
35.27 |
C |
36.53 |
36.05 |
36.46 |
36.20 |
36.92 |
36.41 |
36.87 |
36.42 |
37.14 |
36.57 |
37.12 |
36.56 |
D |
36.69 |
36.26 |
36.63 |
36.34 |
37.22 |
36.61 |
37.18 |
36.64 |
37.39 |
36.80 |
37.36 |
36.65 |
E |
35.97 |
35.12 |
35.93 |
35.54 |
36.50 |
35.47 |
36.41 |
35.73 |
36.82 |
35.79 |
36.82 |
35.90 |
F |
35.56 |
34.72 |
35.50 |
35.08 |
36.11 |
35.08 |
36.03 |
35.26 |
36.34 |
35.28 |
36.30 |
35.47 |
Table 2. SSIM values recovered by different algorithms under different levels of Poisson noise
表2. 不同程度泊松噪声下不同算法复原的SSIM值
|
Peak = 200 |
Peak = 255 |
Peak = 300 |
Image |
Ours |
RLTV |
FOTV |
BM3D |
Ours |
RLTV |
FOTV |
BM3D |
Ours |
RLTV |
FOTV |
BM3D |
A |
0.9677 |
0.9662 |
0.9594 |
0.9617 |
0.9697 |
0.9678 |
0.9604 |
0.9684 |
0.9733 |
0.9718 |
0.9722 |
0.9712 |
B |
0.9569 |
0.9538 |
0.9506 |
0.9521 |
0.9619 |
0.9596 |
0.9610 |
0.9577 |
0.9651 |
0.9645 |
0.9649 |
0.9642 |
C |
0.9674 |
0.9624 |
0.9668 |
0.9650 |
0.9693 |
0.9691 |
0.9691 |
0.9687 |
0.9713 |
0.9725 |
0.9720 |
0.9749 |
D |
0.9640 |
0.9604 |
0.9638 |
0.9502 |
0.9663 |
0.9607 |
0.9666 |
0.9627 |
0.9693 |
0.9614 |
0.9695 |
0.9635 |
E |
0.9599 |
0.9530 |
0.9559 |
0.9486 |
0.9625 |
0.9616 |
0.9612 |
0.9622 |
0.9677 |
0.9662 |
0.9665 |
0.9667 |
F |
0.9599 |
0.9522 |
0.9500 |
0.9402 |
0.9629 |
0.9616 |
0.9621 |
0.9618 |
0.9692 |
0.9667 |
0.9660 |
0.9672 |
从表1和表2中的数据得知,本文的图像复原结果中PSNR与SSIM均有提升,在PSNR数值方面尤为明显。本文算法得到的最大PSNR值比其他算法高出1.37 dB,SSIM值最高高出0.02,在大多数情况下都能达到最佳性能。
本文选取了四张复原图像的细节部分进行展示,对比视觉效果如图1~4所示,四种模型都能很好的复原受损图像。针对图1中Image A图像,我们的模型可以恢复更多纹理的细节,RLTV算法出现锐化过度的问题损害了部分细节信息。针对图2中Image B图像,本文模型在中间黑点处恢复的细节最完整,
Figure 1. Deblurring comparison map of Image A at peak = 255
图1. Image A在噪声peak = 255的去模糊对比图
其他算法都出现大量黑色伪影,边缘细节模糊。针对图3中ImageC图像,BM3D算法复原结果过于平滑,边缘不清晰,不能很好展现图像信息。针对右下角的白点处,FOTV算法丢失了这一细节,本文算法依旧保留了该纹理信息。针对图4中Image D图像,FOTV算法和我们的性能非常相似,没有发现明显的不一样,但是RLTV倾向于放大噪声,噪声没有完全去除干净,残留了噪声信息。BM3D算法过度平滑了表面的纹理信息。由此可见,本文所提出的模型更好地定义了人脑切片图像的边界,在视觉上产生更清晰的结果。
Figure 2. Deblurring comparison map of Image B at peak = 255
图2. Image B在噪声peak = 255的去模糊对比图
Figure 3. Deblurring comparison map of Image C at peak = 200
图3. Image C在噪声peak = 200的去模糊对比图
Figure 4. Deblurring comparison map of Image D at peak = 300
图4. Image D在噪声peak = 300的去模糊对比图
5. 结论
本文设计了小波变换与分数阶全变分模型相结合的医学图像复原算法,对开放的MRI数据集中人脑切片图像复原进行研究。利用小波变换后高频图像具有稀疏性的特征,首先将图像进行小波变换,然后结合分数阶全变分算法,可以有效地去除泊松噪声。使用交替方向乘子法对算法模型进行求解,从而复原出清晰图像。实验表明,我们提出的方法在高精度恢复带有泊松噪声的模糊图像方面具有很大的潜力。未来将考虑加速非盲反卷积算法,提高运算效率。
基金项目
吉林省教育厅科学技术研究项目,JJKH20230788KJ;国家自然科学基金,12171054。