小波域图像降噪

小波域图像去噪

基本模型

考虑加性噪声模型
$y=x+n$y=x+n
其中， y为观测图像， x 为真实图像，n为噪声

等号两边进行小波变换，因为小波变换的叠加性可得
$Y=X+N$Y=X+N
其中 Y 为 y 的小波变换，X 为 x 的小波变换， N 为 n 的小波变换

从统计学的角度来讲，我们要得到真是图像 X 的最有估计 $\hat{X}$X^，即
$\hat{X} = arg \max\limits_X{P(X|Y)}$X^=argXmax​P(X∣Y)
采用贝叶斯模型，即
$P(X|Y)=\frac{P(Y|X)P(X)}{P(Y)}$P(X∣Y)=P(Y)P(Y∣X)P(X)​
则
$\hat{X} = arg \max\limits_X{\frac{P(Y|X)P(X)}{P(Y)}}=arg \min\limits_X \{ -\log P(Y|X)-\log P(X) \}$X^=argXmax​P(Y)P(Y∣X)P(X)​=argXmin​{−logP(Y∣X)−logP(X)}
而
$P(Y|X)=P((X+N)|X)=P(X|X)+P(N|X)=1+P(N)$P(Y∣X)=P((X+N)∣X)=P(X∣X)+P(N∣X)=1+P(N)
因此
$\hat{X} = arg \max\limits_X{\frac{P(Y|X)P(X)}{P(Y)}}=arg \min\limits_X \{ -\log P(N)-\log P(X) \} \tag{1}$X^=argXmax​P(Y)P(Y∣X)P(X)​=argXmin​{−logP(N)−logP(X)}(1)
首先考虑（1）式的第一项

假设噪声 n 服从零均值高斯白噪声，即
$P(n_{ij})\propto exp \{ -\frac{n_{i,j}^2}{2 \sigma_n^2} \}$P(nij​)∝exp{−2σn2​ni,j2​​}
对于正交小波变换，白噪声在小波变换前后性质不变，方差不变，即
$P(N_{i,j})\propto exp \{ -\frac{n_{i,j}^2}{2 \sigma_N^2} \}= exp \{ -\frac{(Y_{i,j}-X_{i,j})^2}{2 \sigma_N^2} \} \tag{2}$P(Ni,j​)∝exp{−2σN2​ni,j2​​}=exp{−2σN2​(Yi,j​−Xi,j​)2​}(2)
在考虑（1）式的第二项，即图像在小波域的统计概率模型，包括广义高斯分布、广义拉普拉斯分布、隐马尔可夫模型等。根据图像的特点，不同的先验统计模型，会得到不同的去噪算法。

这里我们选用最为简单的高斯分布，即
$P(X_{i,j}) \propto exp \{ -\frac{X_{i,j}^2}{2 \sigma_X^2} \} \tag{3}$P(Xi,j​)∝exp{−2σX2​Xi,j2​​}(3)
将（2）、（3）带入到（1）可得：
$\hat{X}_{i,j} = arg \min\limits_{X_{i,j}} \{ \frac{(Y_{i,j}-X_{i,j})^2}{ \sigma_N^2} +\frac{X_{i,j}^2}{ \sigma_X^2} \}$X^i,j​=argXi,j​min​{σN2​(Yi,j​−Xi,j​)2​+σX2​Xi,j2​​}
令
$J= \frac{(Y_{i,j}-X_{i,j})^2}{ \sigma_N^2} +\frac{X_{i,j}^2}{ \sigma_X^2}$J=σN2​(Yi,j​−Xi,j​)2​+σX2​Xi,j2​​
令
$0=\frac{\partial J}{\partial X_{i,j}}=\frac{2（X_{i,j}-Y_{i,j}）}{\sigma_N^2}+\frac{2X_{i,j}}{\sigma_N^2}$0=∂Xi,j​∂J​=σN2​2（Xi,j​−Yi,j​）​+σN2​2Xi,j​​
可得
$\hat{X}_{i,j}= \frac{\sigma_X^2}{\sigma_N^2+\sigma_X^2}Y_{i,j}=WY_{i,j} \tag{4}$X^i,j​=σN2​+σX2​σX2​​Yi,j​=WYi,j​(4)
显然，当选取小波域统计模型为高斯分布时，最优估计的小波域系数为观测图像小波域系数的线性倍数，即维纳滤波。

图像的小波域统计模型为广义高斯分布时的情形

现在考虑小波统计模型为广义高斯分布时的情形

即
$P(X) \propto exp \{ -(\frac{|X|}{\sigma_X})^P\} \tag{5}$P(X)∝exp{−(σX​∣X∣​)P}(5)

P=1\2时

$J_{1/2}= \frac{(Y_{i,j}-X_{i,j})^2}{ \sigma_N^2} +\frac{|X_{i,j}|^{1/2}}{ \sigma_X}\\ 0=\frac{\partial J_1}{\partial X_{i,j}}=0 \Rightarrow \begin{cases} X_{i,j}+\frac{1}{4\sigma_X\sigma_N^2\sqrt{x}}=Y_{i,j} \qquad X_{i,j}\geq 0\\ X_{i,j}+\frac{1}{4\sigma_X\sigma_N^2\sqrt{-x}}=Y_{i,j} \qquad X_{i,j}&lt;0 \end{cases}$J1/2​=σN2​(Yi,j​−Xi,j​)2​+σX​∣Xi,j​∣1/2​0=∂Xi,j​∂J1​​=0⇒{Xi,j​+4σX​σN2​x​1​=Yi,j​Xi,j​≥0Xi,j​+4σX​σN2​−x​1​=Yi,j​Xi,j​<0​

p=1时

$J_1= \frac{(Y_{i,j}-X_{i,j})^2}{ \sigma_N^2} +\frac{|X_{i,j}|}{ \sigma_X}\\ 0=\frac{\partial J_1}{\partial X_{i,j}}=\begin{cases} \frac{2（X_{i,j}-Y_{i,j}）}{\sigma_N^2}+\frac{1}{\sigma_X} \qquad X_{i,j} \geq 0 \\ \frac{2（X_{i,j}-Y_{i,j}）}{\sigma_N^2}-\frac{1}{\sigma_X} \qquad X_{i,j}&lt;0 \end{cases}$J1​=σN2​(Yi,j​−Xi,j​)2​+σX​∣Xi,j​∣​0=∂Xi,j​∂J1​​={σN2​2（Xi,j​−Yi,j​）​+σX​1​Xi,j​≥0σN2​2（Xi,j​−Yi,j​）​−σX​1​Xi,j​<0​

可得
$\hat{X}_{i,j}=\begin{cases} Y_{i,j}-\frac{\sigma_N^2}{2\sigma_X} \qquad X_{i,j} \geq 0\\ Y_{i,j}+\frac{\sigma_N^2}{2\sigma_X} \qquad X_{i,j} &lt;0 \end{cases} \tag{6}$X^i,j​={Yi,j​−2σX​σN2​​Xi,j​≥0Yi,j​+2σX​σN2​​Xi,j​<0​(6)

P=2时

​ 很方便，可以验证最优估计的形式为线性形式，如式（4）。

为了直观地看出取不同P值时，小波域最优估计系数与观测系数的关系，我们做如下图

其中P=1/2时，根据对勾函数的性质可以得到最有估计与观测值的关系如图（a）；P=1时，显然根据式（6）可得图（b）；P=2时，根据式（4）可得图（c）。

事实上，图（a）（b） （c）分别对应了小波域图像去噪收缩模型的硬阈值，软阈值及全局维纳滤波。

在下篇文章中，将针对全局维纳滤波的缺点，介绍局部维纳滤波的方法。