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] = hstar 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的情况这里就不写了.
滤波的推导
自然地, 寻找驻点:
[egin{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(hstar R_{xx}[m] -R_{yx}[m]) \
&=0.
end{array}
]
于是必须满足
[ ag{1}
h star R_{xx}[m] = R_{yx}[m].
]
相应的在频率域内存在(假设DFT存在, 参考文献用的z变换, 这个不是特别了解):
[ ag{2}
H[u] S_{xx}[u] = S_{yx}[u].
]
在FIR情况下, 可以通过(1)推导出一个线性方程组从而求解, non-causal下可用(2)推导出结果.
注: (H[u])用了[]是为了保持一致, 在non-causal情况下(H(z))可能更加妥当.
故最优解为:
[ ag{*}
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)相互独立.
同时
[S_{yy} = S_{yx}.
]
所以:
[ ag{3}
H = frac{S_{yy}}{S_{yy} + S_{vv}}.
]
特别的例子
在图像数字处理中, 给出的这样的情形(FIR):
[x[n] = g star y [n] + v[n],
]
记(r[n] = g star y [n]),
则
[egin{aligned}
S_{rr}[u]
&= mathcal{F}{R_{rr}[m]}[u] \
&= mathcal{F}{E[r[n] r[n-m]]}[u] \
&= E[r[n] mathcal{F}{r[n-m]}][u] \
&= E[r[n] G[-u] Y[-u] e^{-j2pi nu/N}] \
&= G[-u] E[r[n]Y[-u] e^{-j2pi nu/N}] \
&= G[-u] mathcal{F}{R_{ry}[m]}[u] \
&= G[-u] mathcal{F}{E[r[n+m]y[n]]}[u] \
&= G[-u] E[mathcal{F}{r[n+m]}[u]y[n]] \
&= G[-u] E[G[u]Y[u]e^{-j2pi nu/N}y[n]] \
&= G[-u]G[u] mathcal{F}{R_{yy}[m]}[u] \
&= G[-u]G[u] S_{yy}[u]
end{aligned}
]
由于
[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],
]
代替了,
所以
[mathcal{F}{R_{yy}[m]} = F^*(u)F(u) = |F(u)|^2,
]
这里假设(y[n] in mathbb{R}).
(S_{vv}[u])也是类似的.