第9章维纳滤波

第九章维纳滤波

(WienerFiltering)2/5/20231主要内容9.1概述9.2波形线性均方估计的正交原理9.3维纳-霍夫（Wiener-Horf)积分方程9.4非因果的维纳滤波问题9.5因果的维纳滤波器9.6预测问题9.7后验维纳滤波与互补维纳滤波9.8维纳滤波器的应用2/5/20232在实际应用中，有用信号往往会受到一些外界干扰，我们实际观察到的是受到噪声干扰了的信号。如何最大限度地抑制噪声，并将有用信号分离出来，是信号处理中经常遇到的问题。

9.1概述2/5/20233在传输或测量信号s(n)时，由于信道噪声或者测量噪声w(n)，接收或测量到的数据x(n)将与s(n)不同。设噪声是加性的:

即:x(n)=s(n)+w(n)2/5/20234如果s(n)和w(n)的频谱是分离的，那么设计一个具有恰当频率特性的线性滤波器即能有效地抑制噪声并提取信号，这就是前面经典数字信号处理理论中详细讨论过的数字滤波器的设计问题.但是如果s(n)和w(n)的频谱互相重叠,或者s(n)和w(n)是随机信号，它们的频谱根本就不存在，问题就要复杂得多，这就是本章要讨论的内容。2/5/20235随机性是生物医学信号的特点之一，在本章中主要讨论噪声中随机信号的线性估计问题。维纳滤波适用于平稳随机过程。观察x(t)中既含有随机信号s(t)又含有噪声n(t)。2/5/20236经处理器处理后得一估计值作为对所希望取得的信息d(t)的估计值，d(t)可能是s(t)，也可能是预测值s(t+a)，导数ds(t)/dt等。估计的任务就是要求与d(t)的差值在一定判据意义下取得极小值。2/5/20237处理器判据极小9-1估计原理方框图2/5/20238根据待估计量d(t)的形势，波形估计问题可分为三类：滤波问题：由t0~tf一段时期内的观察x(t),t0≤t≤tf，估计t=tf瞬间信号s(t)的值s(t)。即：d(t)=s(t)。预测问题：由t0~tf一段时期内的观察x(t),估计t>tf的某一时刻待估计信号的可能值。即：d(t)=s(t+a)，a>0。平滑问题：由t0~tf一段时期内的观察x(t),估计t0<t<tf期内待估计信号的值.也叫内插问题。2/5/20239在本章应用中作以下假设：x(t)中信号s(t)和噪声信号n(t)是加法地结合的，即：x(t)=s(t)+n(t)。处理器所采用的算法限于线性运算。即：如该处理器的冲击响应是h(t),则：所采用的判据是均方误差最小，即：2/5/2023109.2波形线性均方估计的正交原理设信号是随机过程s(t)，观察是t0~tf期内测得的随机过程x(t),t∊[t0，tf],为了简化讨论，设s(t)和x(t)都是零均值的。采用最小均方误差作为估计判据。即：又限定估计是由观察x(t)经线性滤波器h(t)得出的：2/5/202311最优线性均方估计的选取原则是使估计误差与所有的观察值x(),∊[t0，tf]正交，也就是说，如果对每一个∊[t0，tf]都有：则均方误差最小，它等于：2/5/202312证明：2/5/202313例9.1简单预测问题设观察中没有噪声，即x(t)=s(t)，又待估计量是信号的预测值d(t)=s(t+a)，a>0，设只用t时刻的观察值x(t)对d(t)作线性估计：按最小均方误差判据做估计，即求。2/5/202314解：根据正交原理，估计误差：应和观察值x(t)=s(t)正交，即：2/5/202315例9.2设观察中没有噪声，即x(t)=s(t)，又待估计量是信号的预测值d(t)=s(t+a)，a>0，设估计算子采用：按最小均方误差判据做估计，即求估计系数a和b。2/5/202316解：此时正交原理表现为：2/5/202317由于Rss‘(t)是奇函数，所以Rss‘(0)=0

把上式化简得到：2/5/202318代入均方误差式中，得到：因为Rss(0)>0,此式最后一项大于零，所以，它要比例9.1的LMS要小，主要是他利用了更多的测量信息，估计效果更好些。2/5/2023199.3维纳-霍夫（Wiener-Horf)积分方程观察x(t)由信号s(t)和噪声n(t)相加组成，观察时间t0~tf，则：

x(t)=s(t)+n(t),∊[t0，tf]

待估计过程是d(t)，x(t)经线性处理得到估计为：要求估计均方误差最小，试求h(t)2/5/202320根据正交原理可知：即t时刻的估计误差要和t0~tf区间所有时刻的观察值x()正交，推得：2/5/202321这就是h(t)应满足的条件,称为维纳-霍夫积分方程，只要相关函数Rxd和Rxx已知，就可以由此解出h(t).而h(t)一经解出，就有：问题是维纳-霍夫方程是一个积分方程，未必能求出解析解答。2/5/2023229.4非因果的维纳滤波问题对于滤波问题，利用从t0=-∞一直到tf=t时刻为止的全部观察来估计t时刻的信号s(t)。此时有：

t0=-∞,tf=td(t)=s(t)

于是维纳-霍夫方程变为：2/5/202323做变量替换，t-=,t-=,得到：或：此时：2/5/2023249.4.1连续时间形式的解如果不要求滤波器是因果的，可以把观察时间的上限tf扩展到-∞,也就是利用x(t)在全时间轴上的值来进行估计。此时维纳-霍夫方程为：2/5/202325把它做付氏变换有：如果n(t)和s(t)统计独立，则有：做反傅里叶变换得到滤波器的冲击响应h(t)2/5/202326例9.3

设信号的功率谱是：噪声是白色的，其功率谱是常数而且噪声与信号统计独立，求维纳滤波器的频率特性和冲击响应。

2/5/202327解：此时有求付里叶反变换得到：它显然是非因果的，物理不可实现的th(t)o2/5/202328在离散情况下，在不要求物理可实现的条件下。可以类似推出以下结论：维纳-霍夫原方程为：现在放宽为：9.4.2离散时间形式的解2/5/202329

因此滤波器的频率特性是：实际中一般采用Z变换的传递函数将H(z)做反演Z变换得到冲击响应h(n)2/5/202330可见H(ej)

决定于信号与噪声的功率谱密度;当噪声为零时,即Snn(ej)

=0;H(ej)=1，信号全部通过；当信号为零时,即Sss(ej)

=0;H(ej)=0，噪声被全部抑制掉;因此维纳滤波器确有滤除噪声的能力。非因果维纳滤波器的幅频特性如下图所示。2/5/2023312/5/2023320非因果维纳滤波器的幅频特性wSss(ej)Snn(ej)H(ej)12/5/202333例9.4

设信号的自相关函数是：噪声是白色的设计非因果的维纳滤波器2/5/202334传递函数：解：此时有2/5/202335它是非因果的，而且是无限长的，可以取短近似，如只取4项为：2/5/2023369.5因果的维纳滤波器非因果维纳滤波器需要用全时间上的观察值来估计s(n),所以不能实时实现，即使采用把h(n)截短的近似估计，也必须延迟若干拍，待xn+k输入后（k是截短范围)才能做出本次估计。维纳滤波器的时域解（TimedomainsolutionoftheWienerfilter）2/5/202337设计维纳滤波器的过程就是寻求在最小均方误差下滤波器的单位脉冲响应h(n)或传递函数H(Z)的表达式，其实质就是解维纳－霍夫（Wiener－Horf）方程。我们从时域入手求最小均方误差下的h(n)。这里只讨论因果可实现滤波器的设计，既：

2/5/202338在此主要介绍两种方法：FIR型：限制处理器的形式，只用最近的p+1个观察值[x(n),x(n-1),┅,x(n-p)]来估计s(n),即：预白化处理，把观察序列值x(n)白化。变成白噪声w(n)，再对w(n)做可实现的最优滤波，如图：X(n)d(n)预白化H2w(n)白色2/5/2023399.5.1FIR型处理由正交原理得：2/5/202340令:m=n-m,上式可改写为：用矩阵表示：RxxH=G

2/5/202341只要Rxx是非奇异的，就可以求得H:

H=Rxx-1GRxx是对称且Toplitz型的这时的最小均方误差为：2/5/202342随机信号都可以看成是由一白色噪声W(n)激励一个物理可实现的系统或模型的响应，如图所示：信号模型图9.5.2预白化处理2/5/202343预白化方法是基于如下事实：当x(n)是方差x2=1的白噪声时，有：

Rxx(m-n)=1当m=n=0其它所以上式直接就可以得出：

h(m)=Rxs(m)m=0~∞

离散维纳--霍夫方程为:2/5/202344

对白色的x(n)而言，它的维纳解是Rxs(m)的

m≥0的部分也可以写成:h(m)=Rxs(m)u(m)它的付氏变换写作：

[Sxs(ej)]+或用Z变换表示：

H(z)=[Sxs(z)]+符号[.]+表示原函数m>0部分对应的付氏变换和Z变换。2/5/202345X(n)S(n)+n(n)d(n)=s(n)预白化H2(z)w(n)白色H1(z)E[.]2极小优化-+预白化算法处理框图2/5/202346复习一个Z变换的性质如果h(n)=h(-n),则有:

H(z)=H(z-1);那么，如果z1是Ｈ(Ｚ)的极点，1/Z1一定是H(z)的极点;同样，零点也有这样的性质；还有如果h(n)是实函数，则Ｈ(z)极点一定是共轭对称的．2/5/2023471.预白化滤波器H1(z)的设计

对x(n)的可实现白化滤波器H1(z)可如下求得:2/5/202348

2.最优滤波器H2(z)的设计

因为W(n)是白色的，所以：2/5/202349最后得总滤波器:2/5/202350证明：2/5/202351

例9.5设x(n)=s(n)+n(n),s、n统计独立,且：

设计可实时实现的维纳滤波器。解：2/5/202352因此得：

可见：

2/5/202353括弧中的因子可按下式做部分分式分解：前一项对应与n0部分，后一项对应与n<0部分。分解得A=B=0.6。

2/5/202354因此得:2/5/202355这种处理器有两种具体实现方法：(ⅰ)把h(n)截短成有限长，再与x(n)做卷积：(ⅱ)按IIR方式递归实现：2/5/2023569.6预测问题预测的任务是根据观察x(n)对d(n)的未来值做预报。提法全然不变，只是待估过程d(n)的具体含义变了，因此估计结果也改变。在最常见的一步预测情况下：

因此：2/5/2023579.6.1一般解答因为此时维纳-霍夫方程中的唯一变化是Rxd(m)成为Rxs(m+1),所以一步预测时离散时间的维纳-霍夫方程是:2/5/202358(i)在非因果情况下，把上式放宽成：从而得：当s(n)和n(n)互相独立时，上式变成：2/5/202359（ii）在因果情况下，解答是：例9.6

试为例9.5设计可实现的一步预测器。

因为：2/5/202360本例中：做部分分式分解，上式可分成:前一项对应于n>0部分，故:2/5/202361因此:得:或用IIR方式实现：2/5/2023629.6.2用有限项FIR滤波器实现此时:由正交原理：得维纳-霍夫方程组:2/5/202363用矩阵形式表示:这里的自相关阵仍是对称且Toeplitz型的。2/5/202364例9.7

用p=2的FIR结构给例9.5设计维纳一步预测器。2/5/202365故得:解得：故2/5/2023669.7后验维纳滤波与互补维纳滤波9.7.1后验维纳滤波维纳滤波是以信号和噪声的相关函数或功率谱已有先验知识为前提的。如果这些统计特性未知，就需先作出它们的估计，然后再据以设计维纳滤波器。但是从严格意义上说，这时所得结果并不是真正的维纳滤波器，所以称之为后验“维纳”滤波。2/5/202367从频域上应用后验维纳滤波的核心问题是由各次观察xi(n)中分解出信号的谱估计和噪声的谱估计。通常可采用下述步骤：1.先对各次观察求均值，设做N次观察：式中，s(n)是确定性的诱发响应，ni(n)是第i次刺激后记录中的噪声。2/5/202368

则平均诱发响应：然后求的功率谱。如果s(n)和ni(n)统计独立，各次噪声也互相独立，则:2/5/2023692.再分别对每一次观察xi(n)求功率谱:并求这些功率谱的平均值:2/5/2023703.联立解1,2中最后两式，便可求得Sss(ej)

和Snn(ej)

的估计：2/5/202371据此，得后验"维纳"滤波器如下：（i）用于对单次观察进行滤波：（ii）用于对平均诱发响应x(n)进行滤波：2/5/202372许多研究者用这种滤波方法对各种平均诱发响应进行了滤波，但效果报道不一。有的效果较好，有的却不甚见效。其原因除了谱估计不是真实值，因此所得得H()只能是近似的估计外，还由于“维纳”滤波的其它假设也未必能满足。其中：

过程不是平稳的;“信号和噪声是相加的”这一假设是一个有用模型，但刺激愈接近阈值不正确；信号与噪声未必统计独立。实际上刺激对作为噪声的自发活动往往也有一些作用。2/5/202373为了改进后验“维纳”滤波的效果。又做出了许多改进方案，介绍如下：1.交替集均法此法除按前式求外,又按下式计算另一种平均值：

即：每当序号i为偶数时，就将观察值取负号，通过这样的“相加”，S(n)将被平均掉，因此的功率谱将只反映噪声：2/5/202374便可得"维纳"滤波器。

这个方法的优点式计算量大为下降：只需要求两次功率谱[一次对，一次对]，而采用前法却需要求N+1次功率谱[对每个xi(n)求，还要对求]。但理论分析可以证明，所得谱估计方差较大是这种方法的缺点。2/5/202375

2.谱平滑：把前面得到的功率谱和加以平滑，然后再代维纳滤波公式，可以改善滤波效果。设施加在和上的谱窗口分别是Ｗ1(ej)和Ｗ2(ej)，则平滑后的谱分别是:2/5/202376然后把它们代入维纳滤波公式，得：窗口长度要比较数据长度短，其具体值要在方差和偏差之间取折中。时窗宽则谱窗窄，因此平滑作用小，偏差小，方差大。反之，时窗窄则谱窗宽，平滑作用显著，因此偏差大，方差小。2/5/2023779.7.2互补维纳滤波维纳滤波器得基本假设是信号为随机的，但是实际工作中信号常有些确定性结构，并非纯粹随机，因此应用效果未必好。因此对它简单应用维纳滤波效果未必好，因为这样处理所得的充其量也只是真实S(t)在最小均方误差意义下的逼近，不是真实S(t)。对这类S(t)是确定性信号的情况，采用互补维纳滤波可能是更合理的方案。2/5/202378以做两次观察为例。如果设计滤波器时，H1(z)和H2(z)是分别独立设计的,如图(a)然后再把处理结果相加，效果就未必好.H2(z)H1(z)s(n)+n2(n)s(n)+n1(n)++x(n)2/5/202379如果设计时多引入一个限制条件，如：效果就会好些，如图(b)。1-H1(z)H1(z)s(n)+n2(n)s(n)+n1(n)++x(n)x2(n)x1(n)2/5/202380因为此时:H1(z)s(n)+n2(n)s(n)+n1(n)++x(n)+n1(n)-n2(n)2/5/202381

可见处理结果中信号s(n)成分不变。H1(z)的任务是把n1-n2变成对n2的最优抵消。由于n2和n1-n2都是随机信号,所以应用维纳滤波的效果就比较好。2/5/2023829.8

维纳滤波器的应用（ApplicationofWienerfilter）要设计维纳滤波器必须知道观测信号和估计信号之间的相关函数，即先验知识。如果我们不知道它们之间的相关函数，就必须先对它们的统计特性做估计，然后才能设计出维纳滤波器，这样设计出的滤波器被称为“后验维纳滤波器”。2/5/2023839.8.1.在生物医学信号处理中比较典型的应用就是关于诱发脑电信号的提取。大脑诱发电位（EvokedPotential，EP）指在外界刺激下，从头皮上记录到的特异电位，它反映了外周感觉神经、感觉通路及中枢神经系统中相关结构在特定刺激情况下的状态反应。在神经学研究以及临床诊断、手术监护中有重要意义。2/5/202384EP信号十分微弱，一般都淹没在自发脑电（EEG）之中，从EEG背景中提取诱发电位一直是个难题：EP的幅度比自发脑电低一个数量级，无法从一次观察中直接得到；EP的频谱与自发脑电频谱完全重迭，使得频率滤波失效；在统计上EP是非平稳的、时变的脑诱发电位。通过多次刺激得到的脑电信号进行叠加来提取EP，这是现今