Signals, Systems and Inference, Chapter 11: Wiener Filtering (mit.edu)
基本
在图像处理的时候, 遇到了这个维纳滤波, 其推导的公式不是很理解, 于是上网查了查, 并做个简单的总结.
符号 | 说明 |
---|---|
\(x[k]]\) | 观测信号\(x\)的第k个元素 |
\(\hat{y}\) | 为\(y\)的一个估计 |
\(v\) | 噪声信号 |
\(e[k]\) | 误差, 为\(e[k]=\hat{y}[k] - y[k]\) |
\(R_{xy}(i, j)\) | 相关系数: \(E\{x[n-i]y[n-j]\}\) |
\(S_{xy}\) | 傅里叶变换: \(\mathcal{F} \{R_{xy}\}\) |
基本的假设:
\(y[k]\)服从wide-sense stationary (WSS), 即
- \(E[y[k]] = E[y[0]]\);
- \(R_{yy}(i, j) = R_{yy}[j-i]\).
维纳滤波可以应用于很多场景, 但是这里只讨论下面的去噪的情形:
\[x[n] = y[n] + v[n], \]且假设\(y, v\)之间相互独立, \(E[v]=0\).
我们的目标是找到一个滤波\(h\), 得到一个估计
\[\hat{y}[n] = h\star x [n], \]使得下式最小
\[E[e^2[n]]. \]维纳滤波需要分情况讨论, 这里只关注离散的情形, 包括
- non-causal: \(\hat{y}[n] = \sum_{k=-\infty}^{+\infty} h[k]y[n-k]\);
- FIR: \(\hat{y}[n] = \sum_{k=0}^{N-1} h[k]y[n-k]\).
causal的情况这里就不写了.
滤波的推导
自然地, 寻找驻点:
\[\begin{array}{ll} \frac{\partial E[e^2[n]]}{\partial h[m]} &=2E[e[n]\frac{\partial e[n]}{\partial h[m]}] \\ &=2E[e[n]x[n-m]] \\ &=2(h\star R_{xx}[m] -R_{yx}[m]) \\ &=0. \end{array} \]于是必须满足
\[\tag{1} h \star R_{xx}[m] = R_{yx}[m]. \]相应的在频率域内存在(假设DFT存在, 参考文献用的z变换, 这个不是特别了解):
\[\tag{2} H[u] S_{xx}[u] = S_{yx}[u]. \]在FIR情况下, 可以通过(1)推导出一个线性方程组从而求解, non-causal下可用(2)推导出结果.
注: \(H[u]\)用了[]是为了保持一致, 在non-causal情况下\(H(z)\)可能更加妥当.
故最优解为:
\[\tag{*} H = S_{yx} / S_{xx}. \]特别的情况
进一步, 由于
\[S_{xx} = \mathcal{F}\{R_{xx}\} = \mathcal{F}\{E[(y[n]+v[n])(y[n-m]+v[n-m])]\} =S_{yy} + S_{vv}, \]最后一步成立的原因是\(E{v}=0\), 且\(v, y\)相互独立.
同时
所以:
\[\tag{3} H = \frac{S_{yy}}{S_{yy} + S_{vv}}. \]特别的例子
在图像数字处理中, 给出的这样的情形(FIR):
\[x[n] = g \star y [n] + v[n], \]记\(r[n] = g \star y [n]\),
则
由于
\[S_{yx} = S_{yr} = G[-u]S_{yy}[u], \]证明和上面是类似的,
\[S_{xx} = S_{rr} + S_{vv}, \]于是
\[H[u] = \frac{1}{G[u]}\frac{1}{1 + S_{vv}[u] / (G[-u]G[u]S_{yy}[u])}. \]当\(g\)为实的时候, 有\(G[-u]G[u]=|G[u]|^2\), 在数字图像处理书中, 给出的公式中:
\[S_{yy}[u] = |F[u]|^2, S_{vv}[u] = |V[u]|^2, \]个人觉得是这里的期望\(E\)用
\[R_{yy}[m] = \sum_{n=0}^{N-1} y[n]y[n-m], \]代替了,
所以
这里假设\(y[n] \in \mathbb{R}\).
\(S_{vv}[u]\)也是类似的.