描述
开 本: 16开纸 张: 胶版纸包 装: 平装-胶订是否套装: 否国际标准书号ISBN: 9787302492351丛书名: 高等学校电子信息类专业系列教材
本书力求以简约的内容概括随机信号处理的基本理论与方法,为读者打下牢固的随机信号处理知识基础。本书既可作为电子信息类高年级本科生和相关学科研究生的教材,也可为从事相关领域研究的科研人员提供参考。
目录
第1章数字信号处理基本概念
1.1概述
1.2离散时间信号
1.2.1连续时间信号的采样
1.2.2采样定理
1.2.3几种常见的数字信号
1.2.4信号的能量、功率及周期性
1.2.5信号的基本运算
1.3信号的傅里叶变换
1.3.1连续时间信号的傅里叶变换
1.3.2离散时间信号的傅里叶变换
1.3.3离散傅里叶变换及其性质
1.4z变换
1.4.1z变换的定义
1.4.2z变换的收敛域
1.4.3z变换的性质
1.4.4逆z变换
1.5离散时间系统
1.5.1基本概念
1.5.2离散时间系统的单位冲激响应函数
1.5.3LSI系统的稳定性和因果性
1.5.4LSI系统的变换域分析
本章习题
第2章随机信号分析基础
2.1概述
2.1.1随机信号的基本概念
2.1.2随机信号的分类
2.2随机信号的概率结构
2.2.1概率论基本概念
2.2.2随机信号有限维概率密度及数字特征
2.3随机信号的平稳性
2.4离散随机信号和复随机信号
2.4.1离散时间随机信号及其数字特征
2.4.2复随机信号
2.5随机信号的遍历性
2.5.1总集意义上的数字特征与时间意义上的数字特征
2.5.2平稳随机信号的遍历性
2.6平稳随机信号的功率谱密度
2.6.1维纳辛钦定理
2.6.2功率谱密度的性质
2.6.3离散随机序列的功率谱密度
2.7几种常见的随机信号
2.7.1白噪声
2.7.2高斯随机信号
2.7.3马尔可夫随机信号
本章习题
第3章信号估计与检测基础
3.1估计的基本概念
3.2估计算法的性能指标
3.2.1性能指标
3.2.2随机信号均值及相关函数的估计
3.3估计性能界——CRB
3.3.1单参数实常量估计的CRB
3.3.2多参量估计的CRB
3.3.3参数变换的CRB
3.3.4复参数估计的CRB
3.4最大似然估计
3.4.1最大似然估计的基本原理
3.4.2变换参数的最大似然估计
3.5贝叶斯估计
3.5.1代价函数
3.5.2最小均方误差估计
3.5.3条件中位数估计
3.5.4最大后验概率估计
3.5.5贝叶斯估计举例
3.6线性最小均方误差估计
3.6.1随机参量的线性最小均方误差估计
3.6.2线性最小均方误差估计的几何解释
3.7最小二乘估计
3.8信号检测基础
3.8.1确定性信号检测
3.8.2随机信号检测
本章习题
第4章随机信号的更新与建模
4.1随机信号通过线性系统
4.1.1基本概念
4.1.2线性系统输入输出信号之间数字特征的关系
4.2随机矢量的线性变换
4.3离散时间序列的线性模型
4.3.1离散时间序列的自回归滑动平均模型
4.3.2ARMA模型的传递函数
4.3.3ARMA系统的等效性
4.4ARMA模型的数字特征
4.4.1互相关函数
4.4.2自相关函数
4.4.3功率谱密度
4.5ARMA、AR、MA模型之间的关系
4.5.1Wold分解定理
4.5.2柯尔莫可洛夫定理
本章习题
第5章随机信号的滤波
5.1数字滤波器的基本概念
5.2维纳滤波
5.2.1最小均方误差(MMSE)准则与正交性原理
5.2.2WienerHopf正则方程
5.2.3Wiener滤波器的求解
5.3线性预测
5.4卡尔曼滤波
5.5最小二乘滤波
5.6匹配滤波器
5.7自适应滤波
5.7.1自适应滤波器的基本概念
5.7.2LMS自适应滤波器
5.7.3RLS自适应滤波
本章习题
第6章功率谱估计
6.1概述
6.2经典谱估计的基本方法
6.2.1经典谱估计法一——周期图法
6.2.2经典谱估计法二——间接法(相关图法、BT谱估计)
6.2.3经典谱估计方法的改进
6.3功率谱估计的参数模型法
6.3.1AR谱估计的相关函数法
6.3.2LevinsonDurbin算法
6.3.3AR谱估计的性质
6.3.4MA谱估计、ARMA谱估计
6.4特征分解法谱估计
6.4.1Pisarenko谐波分解与相关矩阵的特征分解
6.4.2子空间法功率谱估计
本章习题
前言随机信号处理涉及多媒体信号、生物医学信号、通信信号、控制信号等方面,有非常广阔的研究范围。本书力图以有限的篇幅对随机信号处理的基本理论和方法进行概括,使学生对随机信号概念和数学建模有必要的了解,掌握随机信号理论和分析处理的基本方法。全书参考学时为36学时,如需补充相关学科的前沿知识,可拓展到48学时。本书需要“概率论”“信号与系统”“数字信号处理”课程作为预备知识。在内容编排上,本书兼顾实信号和复信号,所以还需要一定的复变函数知识。作者在编写本教材的过程中,参考了大量的国内外学者文献,在此对相关作者表示衷心的感谢。同时,本书的编写也得到清华大学出版社的大力支持,在此表示诚挚的感谢。由于作者水平有限,本书在内容的选择、体系的安排以及文字叙述上难免有疏漏,恳请读者批评指正。编者2018年4月
信号估计与检测基础
3.1估计的基本概念信号处理的很多问题以参数估计的形式存在,例如在雷达及声呐信号处理中,需要估计目标(飞机或船只)的距离和到达角度,这些参数可以通过对发射信号和回波信号的处理获得。又如在现代的高保真录制系统中,需要精确保持录音基准电平,这一般通过电平自动控制系统实现,具体流程是周期性地估计该基准电平,并根据估计结果做出调整。信号处理中的参数估计一般基于观测信号与待估计参数之间的依赖关系,在不同的问题模型中,这种依赖关系可能是确定性的或随机性的。例如在通信中的信道估计问题中,发送信号x(n)经过冲激响应为h(n)的信道输出为y(n)=h(n)x(n)(31)
如果系统采用基于训练序列的信道估计方法,x(n)为已知,这时信道输出y(n)与信道参数h(n)存在确定性依赖关系; 如果系统信道估计不采用训练序列,而是基于信道输出的统计特性,则y(n)与信道参数h(n)存在随机性依赖关系,具体而言,就是y(n)的概率分布函数与h(n)存在依赖关系。随机性依赖关系与确定性依赖关系相比,更贴近实际环境。这首先是因为实际环境中观测信号的获取总是受到干扰,如环境干扰以及观测设备自身的热噪声干扰等。考虑噪声影响后的信道接收信号一般表示为y(n)=h(n)x(n) v(n)(32)其中v(n)表示系统环境干扰或设备热噪声产生的加性干扰,一般情况下用零均值、独立同分布的高斯随机信号描述。这时y(n)也体现为具有高斯分布特性的随机信号,h(n)对y(n)的影响也体现为对其概率密度分布和统计特性的影响。随机性依赖关系更具重要性的另一原因是实际系统中待估计参数经常具有时变性,例如雷达信号处理中目标距离就是不断变化的参数; 移动通信中,发送端或接收端的移动会造成信道环境的改变。待估计参数的变化不是毫无规律,而是服从一定的概率密度分布,例如非视距移动信道的时变性服从瑞利分布。在式(31)中,如果h(n)是具有一定概率结构的随机信号,则无论x(n)是确定性信号或者随机信号,h(n)对y(n)的影响都体现为概率密度分布和统计特性的影响。下面以一个简单的例子来说明参数估计的过程。例31在电平自动控制系统的电平估计问题中,对某个电平值θ进行估计,由于环境干扰及测量设备的不完善,测量总会存在误差,测量误差可归结为噪声,因此,实际得到的测量值为x=θ w(33)其中,w一般服从零均值高斯分布,方差为σ2。问题是,如何根据测量值x来估计θ的值。一种最简单的方法是将多个观测结果进行平均。首先定义一个时间窗口,认为实际电平值在这个时间窗口中是恒定的。在这个时间窗口里得到电平值的N个含噪测量结果,记为x(0),x(1),…,x(N-1)。这些测量结果中x(n)=θ w(n),n=0,1,…,N-1(34)
根据x(n)得到θ的估计为θ^=x(0) x(1) … x(N-1)N(35)
上式的一个直观效果是通过对噪声求均值,以降低噪声的影响。一般情况下,可以认为估计结果与观测数据由一个映射关系确定: θ^=f(x(0),x(1),…,x(N-1))(36)
上面提到,由于噪声的影响,x(n)更多时候体现为随机信号的特征,因此上式中x(0),x(1),…,x(N-1)可看成随机变量。θ^与x(0),x(1),…,x(N-1)存在映射关系,可认为θ^也是随机变量(不管θ是常量或随机变量)。一次估计的结果可以看成是θ^的一个样本,因此,验证估计算法的性能往往需基于多次估计结果进行评估。实际应用中往往很难得到类似于式(35)的闭式表达式,参数估计可能包含复杂的步骤,但包含的要素是一致的,从上面这个例子可以看出构造一个估计问题的基本要素如下: (1) 待估计参数: 指估计的目标参数,可以是单参量,也可以是多参量。待估计的参数可能是常量,也可能是随机变量,如果是随机变量一般要求其概率密度函数已知。(2) 信号模型: 指观测信号的产生模型,其中必须体现待估计参数与观测信号的依赖关系。观测信号可以是单参量,也可以是多参量。(3) 误差准则: 指估计算法的性能衡量指标,很多时候误差准则可指导估计算法的设计。误差准则的确定也有助于不同估计算法间的比较,特别是针对同一个估计问题的不同估计算法。(4) 估计算法: 对同一估计问题,估计算法不是唯一的,不同的估计算法可以得到不同的估计性能,如何得到误差小、计算量低的估计算法是信号处理的核心问题。3.2估计算法的性能指标依据不同的估计准则,可以得到不同的估计量,这些估计量的性能需要有评价的指标。估计量是观测的函数,而观测是随机变量(或矢量),因此,估计量也是随机变量,随着观测数据的增多,希望估计量能逐渐逼近真值,因此,对估计量的均值和方差应有一定的要求。3.2.1性能指标估计量的好坏可以从无偏性、有效性、一致性来加以评价。1. 无偏性待估计量为常数θ,它的估计值为θ^,如果E(θ^)=θ,则称θ^为θ的无偏估计,否则称θ^为有偏估计。定义估计的偏差为θ~=E(θ^)-θ(37)
通常希望估计量的均值趋于被估计量的真值或被估计量的均值,即估计应该是无偏的。如果估计θ^不是无偏估计,但随着样本数目的增加,其数学期望趋近于真的估计量,即limN→∞[E(θ^)-θ]=0(38)则称估计θ^是渐近无偏估计。当被估计量θ是随机变量时,如果估计量的均值等于被估计量的均值,即E[θ^]=E[θ](39)
则称θ^为无偏估计。对于有偏估计,如果观测数据越多,估计的性能越好,即limN→∞[E(θ^)-E(θ)]=0(310)则称θ^为随机变量的渐近无偏估计。2. 有效性估计量具有无偏性并不表明已经保证了估计的品质。当被估计量为未知常数时,不仅希望估计量的均值等于真值,而且希望估计量的取值集中在真值附近,这一品质可以通过估计的方差来描述,估计的方差为Var(θ^)=E{[θ^-E(θ^)]2}(311)对于无偏估计,方差越小,表明估计量的取值越集中,估计的性能越好,估计也越有效。对于有偏估计,估计的方差小并不能说明估计是好的,因为如果估计有偏差,方差小的估计仍然可能有较大的估计误差,这时用均方误差加以描述更合适,估计的均方误差定义为Mse(θ^)=E{|θ^-θ|2}(312)均方误差越小,表明估计越有效。3. 一致性当用N个观测值估计参量时,一般来说,观测值越多,估计越趋于真值,如果limN→∞P{θ-θ^=1N2∑N-1n=0∑N-1m=0E[x(n)x(m)]
=1N2∑N-1n=0E[x2(n)] ∑N-1n=0∑N-1m=0,m≠nE[x(n)x(m)]
=1N2∑N-1n=0E[x2(n)] ∑N-1n=0∑N-1m=0,m≠nE[x(n)]E[x(m)](316)
记E[x2(n)]=Dx,E[x(n)]=mx,所以E[m^2x]=1NDx N(N-1)N2m2x
=1NDx N-1Nm2x(317)所以,估计的方差为Var(m^x)=E(m^2x)-[E(m^x)]2
=1NDx N-1Nm2x-m2x=1Nσ2x(318)可见,当N→∞时,估计的方差趋于零,所以m^x是一致估计。2. 自相关函数的均值在实际工作中,对平稳随机序列x(n),所能得到的x(n)的N个观察值为x(0),x(1),…,x(N-1),记为xN(n)。对n≥N时的x(n)的值假设为零。现在的任务是如何由这N个观察值来估计出x(n)的自相关函数Rx(m)。对Rx(m)的估计方法通常有两种,一是直接计算法,二是利用FFT来求解。下面以自相关函数的直接估计法为例分析其性能。如果观察值的点数N为有限值,则求R(m)估计值的直接估计法是R^(m)=1N∑N-1n=0xN(n)xN(n-m)(319)
由于x(n)只有N个观察值,因此,对于每一个固定的迟延m,在0~N-1的范围内可以利用的数据只有N-|m|个,所以在实际计算R^(m)时,式(319)变为R^(m)=1N∑N-1-|m|n=0x(n)x(n-m)(320)
值得注意的是,R^(m)的长度为2N-1,对实信号,它是以m为对称中心的偶对称函数。所以在上式的计算中,可先计算m≤0的部分,再根据对称性求m>0的部分。现在讨论R^(m)估计的质量。1) 估计的偏差由偏差的定义有Bia[R^(m)]=E[R^(m)]-R(m)(321)
式中E[R^(m)]=E1N∑N-1-|m|n=0x(n)x(n-m)
=1N∑N-1-|m|n=0E[x(n)x(n-m)]
=1N∑N-1-|m|n=0R(m)即E[R^(m)]=N-|m|NR(m)(322)
所以估计的偏差为Bia[R^(m)]=-mNR(m)(323)
分析式(322)和式(323),可以看出,估计不是无偏估计,而且当m→N时,估计的偏差很大。对有限的m,当N→∞时,有E[R^x(m)]=Rx(m),所以是渐进无偏的。由此,可以得到如下结论: (1) 对于一个固定的延迟m,当N→∞时,Bia[R^(m)]→0,因此,R^(m)是对R(m)的渐近无偏估计; (2) 对于一个固定的N,只有当|m| <0,|m|≥N(324)
图31三角窗w(m)
的乘积,w(m)的长度是2N-1,如图31所示。式(324)的三角窗函数又称Bartlett窗,由于它对R(m)进行了加权,致使R^(m)产生了偏差。显然,由于加权窗的存在,产生了上述第(2)点结论。
该窗函数实际上是由于对真实数据的截短而产生的。因为xN(n)可以看作真实数据x(n)和一矩形窗函数d(n)相乘的结果,即xN(n)=x(n)d(n)(325)式中d(n)=1,0≤n≤N-1
0,n<0,n≥N(326)
根据式(320)和式(325)可得R^(m)=1N∑N-1-|m|n=0xN(n)xN(n-m)
=1N∑N-1-|m|n=0x(n)d(n)x(n-m)d(n-m)(327)E{R^(m)}=1N∑N-1-|m|n=0E{x(n)x(n-m)}d(n)d(n-m)
=R(m)N∑N-1-|m|n=0d(n)d(n-m)
=R(m)w(m)(328)
可见,w(m)正是对矩形数据窗函数d(n)作自相关的结果。由上面的讨论可知,当对一个信号做自然截短时,就不可避免地对该数据施加了一个矩形窗口,由此数据窗口就产生了加在自相关函数上的三角窗口。三角窗口影响了R^(m)对R(m)估计的质量。加在数据上的窗口,一般称为数据窗。加在自相关函数上的窗口一般称为延迟窗,后面将看到,这些窗函数也直接影响了功率谱估计的质量。2) 估计的方差根据方差的定义,有Var[R^(m)]=E[R^(m)-E[R^(m)]]2
=E[R^2(m)]-{E[R^(m)]}2(329)由式(322)可得{E[R^(m)]}2=N-|m|NR(m)2(330)同时可得E[R^2(m)]=E1N2∑N-1-|m|n=0x(n)x(n-m)∑N-1-|m|k=0x(k)x(k-m)
=1N2∑n∑kE[x(n)x(k)x(n-m)x(k-m)](331)
这里要计算随机信号x(n)的四阶矩,对一般平稳信号的四阶矩很难计算,假定x(n)是零均值的高斯随机信号,可以证明,其四阶矩满足如下关系式E[x(n)x(k)x(n-m)x(k-m)]
=E[x(n)x(k)]E[x(n-m)x(k-m)]
E[x(n)x(n-m)]E[x(k)x(k-m)]
E[x(n)x(k-m)]E[x(k)x(n-m)]
=R2(n-k) R2(m)
R(n-k m)R(k-n m)(332)
所以E[R^2(m)]=1N2∑n∑k[R2(n-k) R2(m) R(n-k m)R(k-n m)]
=N-|m|NR(m)2 1N2∑n∑k[R2(n-k)
R(n-k m)R(k-n m)](333)将式(330)和式(333)代入式(329),可得Var[R^(m)]=1N2∑N-1-|m|n=0∑N-1-|m|k=0[R2(n-k) R(n-k m)R(k-n m)](334)对任意序列g(n),下式成立∑N-1-|m|n=0∑N-1-|m|k=0g(n-k)=∑N-1-|m|i=-(N-1-|m|)(N-|m|-|i|)g(i)(335)
因此令n-k=i,可把式(334)的双求和变成单求和,即对于一个确定的i,n-k=i的组合数为N-|m|-|i|。
Var[R^(m)]=1N∑N-1-|m|i=-(N-1-|m|)1-|m| |i|N[R2(i) R(i m)R(i-m)](336)
显然,当N→∞时,Var[R^(m)]→0,又因为limN→∞Bia[R^(m)]=0,所以,可得如下结论: 对固定的延迟m,R^(m)是R(m)的一致估计。对R(m)的另一种直接估计方法是对式(320)稍做修改,即R^(m)=1N-|m|∑N-1-|m|n=0xN(n)xN(n m)(337)可以证明,该式的R^(m)是对R(m)的无偏估计,利用上述类似的推导,可分析无偏估计子的性能,最终可得有偏估计子: Var[R^(m)]∝O1N
无偏估计子: Var[R^(m)]∝O1N-|m|
可以看到无偏估计子的方差性能较差。实际应用一般采用有偏估计子,除了方差的原因外,还基于以下考虑: ①当m取值较小时,有偏估计子的偏差不明显; ②对大多数平稳随机信号,当m值较大时,R(m)的值一般很小。式(323)显示偏差的值也和R(m)的取值有关系,而m较大时R(m)的取值一般会变小,这时可以认为偏差也比较小。3.3估计性能界——CRB不同的估计方法可以得到不同的估计量,估计量的性能可以通过前面介绍的无偏性、有效性和一致性来评价,但在实际中,估计量可能比较复杂,很难评价估计量的有效性和一致性。此外,在得到一个估计量后,它的性能是否已经达到最佳?是否还有更好的估计量?克拉美罗下限(CramérRao Bound,CRB)揭示了无偏估计量估计方差的最小值。3.3.1单参数实常量估计的CRB首先考虑只有一个待估计参数的情况,令θ为待估计参数,并且为实参数(非随机量)。观测信号x(n)与θ的依赖关系用条件概率p(x|θ)表示(为描述方便,以下将x(n)简记为xn),有如下定理。定理31令x=(x0,x1,…,xN-1)为观测样本矢量,p(x|θ)是x的条件概率密度。若θ^是θ的一个无偏估计子,且f(x|θ)/θ存在,则θ^的方差满足下式Var(θ^)=E[(θ^-θ)2]≥I-1(θ)(338)
式中I(θ)=Elnp(x|θ)θ2=-E2lnp(x|θ)θ2(339)当且仅当lnp(x|θ)θ=K(θ)(θ^-θ)(340)
时,式(338)的等号成立。其中K(θ)是θ的正函数,且不包含x。证明: 由假设条件知E(θ^-θ)=0,因此E[θ^-θ]=∫ ∞-∞…∫ ∞-∞(θ^-θ)f(x|θ)dx0…dxN-1
=Δ∫ ∞-∞(θ^-θ)f(x|θ)dx
=0(341)
上式两边对θ求偏微分可得θE[θ^-θ]=∫ ∞-∞θ[(θ^-θ)f(x|θ)]dx=0(342)由上式得∫ ∞-∞-f(x|θ)dx ∫ ∞-∞(θ^-θ)θf(x|θ)dx=0(343)
因为p(x|θ)是概率密度函数,所以∫ ∞-∞f(x|θ)dx=1(344)
同时根据组合函数导数性质可得θf(x|θ)=θlnf(x|θ)·f(x|θ)(345)
式(343)重写为∫ ∞-∞θlnf(x|θ)·(θ^-θ)f(x|θ)dx=1(346)
上式进一步重写为∫ ∞-∞θlnf(x|θ)f(x|θ)·(θ^-θ)f(x|θ)dx=1(347)根据泛函分析中的柯西不等式|〈z1·z2〉|2≤〈z1·z1〉·〈z2·z2〉,其中z1、z2是泛函空间的元素,〈·〉表示内积操作。式(347)中令z1=f1(x)=lnf(x|θ)/θf(x|θ)z2=f2(x)=(θ^-θ)f(x|θ)
内积操作定义为∫ ∞-∞f1(x)f2(x)dx。由式(347)得∫ ∞-∞θlnf(x|θ)2f(x|θ)dx·∫ ∞-∞(θ^-θ)2f(x|θ)dx≥1(348)等价于∫ ∞-∞(θ^-θ)2f(x|θ)dx≥1∫ ∞-∞θlnf(x|θ)2f(x|θ)dx(349)上式左边刚好就是无偏估计的方差计算结果。柯西不等式等于号成立的条件是z1,z2成比例,即式(340)成立。证毕特别要指出的是式(342)积分和求导顺序可以互换,这在实际中一般可以满足。为了得到更直观的描述,将式(344)两边对θ求导数,并将积分和求导顺序互换可得∫ ∞-∞θf(x|θ)dx=0(350)
基于式(345)可得∫ ∞-∞θlnf(x|θ)·f(x|θ)dx=Eθlnf(x|θ)=0(351)
上式称为定理3.1成立的规范化条件,是积分和求导顺序可以互换的直观描述。该规范化条件两边对θ求导数,则可得式(339)。CRB给出了无偏估计量估计方差的下限。达到CRB的估计,其估计的方差是最小的,称这样的估计为有效估计。I(θ)称为数据的费希尔(Fisher)信息,CRB是费希尔信息的倒数。直观理解是,费希尔信息越多,CRB越低。费希尔信息具有信息度量的基本性质。首先,由式(339)可以看出,费希尔信息具有非负性; 其次,对于独立观测,费希尔信息具有可加性,这是因为,对于独立的观测有lnp(x|θ)=∑N-1n=0lnp(xn|θ)(352)
所以-E2lnp(x|θ)θ2=-∑N-1n=0E2lnp(xn|θ)θ2(353)如果观测是独立同分布的,那么I(θ)=Ni(θ)(354)
式中i(θ)=-E2lnp(xn|θ)θ2(355)
对于非独立观测,费希尔信息要比独立观测低,例如,对于x0=x1=…=xN-1这种完全相关的情况,I(θ)=i(θ),此时,CRB并不随观测数据的增加而降低。例32独立同分布平稳随机信号服从高斯分布,用式(314)估计该随机信号的均值,该估计子的方差是否达到CRB?解平稳高斯信号的一维概率分布为p(x|mx)=12πσexp-(x-mx)22σ2(356)多变量的联合概率密度为p(x|mx)=12πσN∏N-1n=0exp-(xn-mx)22σ2(357)对数联合概率密度的一、二阶偏导数为lnp(x|mx)mx=∑N-1n=0xn-mxσ2(358)2lnp(x|mx)m2x=-Nσ2(359)可得到其费希尔信息量为I(mx)=-E2lnp(x|mx)m2x=Nσ2(360)
对比3.2.2节中估计子的方差,可以看到该均值估计可达到CRB。注意到式(358)中的一阶导数可重写为lnp(x|mx)mx=Nσ2(m^x-mx)(361)即一阶导数和(m^x-mx)成正比,显示该估计子满足式(340)的条件。例33信号观测模型为xn=Acos(2πf0n ) w(n),其中f0为已知参数,w(n)∝N(0,σ2),相位为未知参数。求相位估计的CRB。解观测信号的一维概率分布为p(x|)=12πσexp-(xn-Acos(2πf0n ))22σ2(362)假设参数估计基于N个观测,记为x=[x0, x1, …, xN-1]T,其联合概率函数为p(x|)=p(x0,x1,…,xN-1|)
=12πσN∏N-1n=0exp-(xn-Acos(2πf0n ))22σ2(363)
可得到观测信号的对数似然函数为lnp(x|)=-N2ln(2πσ2)-12σ2∑N-1n=0(xn-Acos(2πf0n ))2(364)
可得lnp(x|)=-1σ2∑N-1n=0[xn-Acos(2πf0n )]Asin(2πf0n )
=-Aσ2∑N-1n=0[xnsin(2πf0n )-A2sin(4πf0n 2)]2lnp(x|)2=-Aσ2∑N-1n=0[xncos(2πf0n )-Acos(4πf0n 2)]-E[2lnp(x|)2]=Aσ2E∑N-1n=0[xncos(2πf0n )-Acos(4πf0n 2)]
=Aσ2∑N-1n=0E[xn]cos(2πf0n )-Acos(4πf0n 2)
=Aσ2∑N-1n=0Acos2(2πf0n )-Acos(4πf0n 2)
=A2σ2∑N-1n=012 12cos(4πf0n 2)-cos(4πf0n 2)
=A22σ2∑N-1n=01-cos(4πf0n 2)=NA22σ2-A22σ2∑N-1n=0cos(4πf0n 2)(365)
上式中求和项为周期序列,有limN→∞1N∑N-1n=0cos(4πf0n 2)=0(366)当N比较大时,可认为1N∑N-1n=0cos(4πf0n 2)<<1(367)
于是可认为式(365)中第二项远小于第一项,直接忽略第二项可得-E2lnp(x|)2=NA22σ2(368)其倒数即为相位估计的CRB。3.3.2多参量估计的CRB假设有L个未知参数,用矢量表示为θ=[θ1,θ2,…,θL]T(369)
参数估计的CramerRao不等式表示为E[(θ^-θ)(θ^-θ)T]≥I-1(θ)(370)
其中符号“≥”用于矩阵时有不同的含义。对矩阵A、B,A≥B表示A-B是半正定矩阵,A>B表示A-B是正定矩阵。上式中左边为估计子的协方差矩阵,右边的I(θ)称为费希尔信息矩阵,定义为I(θ)=-E2lnp(x|θ)θθT(371)
上式中用到了矢量函数的求导公式,下面给出其定义。设f(θ)是L×1矢量θ的函数,函数值为标量,其导数表示为f(θ)θ=f(θ)θ1f(θ)θ2 f(θ)θLL×1(372)
根据该定义可得以下常用导数aTθθ=θTaθ=a(373)θTAθθ=(AT A)θ(374)式(372)的定义可拓展到多个函数的组合,假设有R个关于θ的函数构成如下矢量: f(θ)=f1(θ)f2(θ)fR(θ)R×1(375)则其关于θ的导数定义为f(θ)θT=f1(θ)θ1f1(θ)θ2…f1(θ)θL
f2(θ)θ1f2(θ)θ2…f2(θ)θL
fR(θ)θ1fR(θ)θ2…fR(θ)θLR×L(376)
根据上述定义,可得矩阵I(θ)中第{i,j}个元素为[I(θ)]ij=-E2lnp(x|θ)θiθj(377)对式(371)这里不做详细证明,下面以一个例子介绍多参量CRB的计算。例34对独立同分布平稳高斯随机信号,计算均值估计和方差估计的CRB。解由于待估计矢量为θ=[mxσ2]T,因此,费希尔信息矩阵是2×2的矩阵,即I(θ)=-E2lnp(x|θ)θθT=-E2lnp(x|θ)m2x-E2lnp(x|θ)mxσ2
-E2lnp(x|θ)σ2mx-E2lnp(x|θ)(σ2)2(378)
由式(356)得到对数似然函数为lnp(x|θ)=-N2ln(2πσ2)-12σ2∑N-1n=0(xn-mx)2(379)
所以lnp(x|θ)mx=1σ2∑N-1n=0(xn-mx)(380)lnp(x|θ)σ2=-N2σ2 12σ4∑N-1n=0(xn-mx)2(381)2lnp(x|θ)m2x=-Nσ2(382)2lnp(x|θ)(σ2)2=N2σ4-1σ6∑N-1n=0(xn-mx)2(383)2lnp(x|θ)mxσ2=-1σ4∑N-1n=0(xn-mx)(384)2lnp(x|θ)σ2mx=2lnp(x|θ)mxσ2(385)
费希尔信息矩阵为I(θ)=-E-Nσ2-E-1σ4∑N-1n=0(xn-mx)
-E-1σ4∑N-1n=0(xn-mx)-EN2σ4-1σ6∑N-1n=0(xn-mx)2=Nσ200N2σ4(386)
因此,容易得出参数估计的CRB为Var(m^x)≥σ2N(387)
Var(σ^2)≥2σ4N(388)
3.3.3参数变换的CRB在实际中,经常遇到希望估计的参数是某个参数的函数,例如,在电平估计中,有时感兴趣的不是信号电平A的估计,而是信号的功率A2,假设估计A的CRB已知,如何得到A2估计的CRB呢?以下公式回答这一问题。假定可以直接计算CRB的参数为θ,希望估计的参数为α,两者有如下函数关系: α=T(θ)。那么,可以证明,α^的CRB为Var(α^)≥T(θ)θ2-E2lnp(x|θ)θ2(389)
例如,基于例33中式(387)的结果求α=m2x的CRB,则Var(α^)≥(2mx)2N/σ2=4m2xσ2N(390)实际应用中经常碰到多变量的映射关系。例如,独立同分布高斯随机信号的平均功率D=m2x σ2。前面已经得到估计mx、σ2的CRB,在此基础上可采用以下方法得到估计D的CRB。假设P×1维矢量θ包含模型中的可估计参数,有多个待估计参数,并且均可表示为θ的函数,记为α=T(θ),其中T(θ)是R×1维的函数矢量,即α=[T1(θ)T2(θ)…TR(θ)]T(391)
则E[(α^-α)(α^-α)T]≥T(θ)θTI-1(θ)T(θ)θTT(392)式中T(θ)/θT的定义见式(376)。例35对高斯白噪声中恒定电平,求信噪比估计的CRB。解高斯白噪声中恒定电平观测的数学模型为xn=A wn,其中A为电平值,wn∝(0,σ2)为零均值高斯白噪声,信噪比定义为α=A2/σ2,因此可先计算估计A,σ2的CRB。注意到观测数据可以看成是非零均值高斯信号,A,σ2的估计对应均值和方差的估计,可采用例34的结果。设θ=[Aσ2]T,T(θ)=A2/σ2,由例34可知I(θ)=Nσ200N2σ4(393)又因为T(θ)θT=T(θ)AT(θ)σ2=2Aσ2-1σ4(394)
所以T(θ)θTI-1(θ)T(θ)θTT=2Aσ2-1σ4σ2N002σ4N2Aσ2-1σ4(395)=4A2Nσ2 2Nσ4(396)
可以证明,对矢量参数的有效估计量经过线性变换后仍然是有效估计量,而矢量参数的有效估计量经过非线性变换后是渐近有效估计量。3.3.4复参数估计的CRB当待估计参数是复数时,推导涉及复变函数的微分操作,比实参数估计复杂很多。这一节先给出复变函数的基本概念,再介绍一种复参数估计CRB的推导方法。令复自变量为z=x jy,x,y∈R。定义在复数域的函数为f(z)=u(z,z) jv(z,z)=u(x,y) jv(x,y)(397)
其中u(x,y),v(x,y)为实二元函数。f(z)的导数定义为fz=12fx-jfy
fz=12fx jfy(398)复变函数可导的充分必要条件可参考相关复变函数文献。这里假设所涉及的复变函数均可导,以下是一些常用的性质。假设g(z),f(z)在区域D内可导,则其导数满足以下性质: ddz(f(z) g(z))=df(z)dz dg(z)dz(399)ddz(f(z)·g(z))=df(z)dzg(z) dg(z)dzf(z)(3100)ddzf(z)g(z)=1g2(z)df(z)dzg(z)-dg(z)dzf(z)(3101)ddz(f(g(z)))=df(g)dgdg(z)dz(3102)
对多变量函数,如f(z),其中z=[z1,z2,…,zN]T为包含所有自变量的矢量,其导数定义为f(z)z=f(z)z1,f(z)z2,…,f(z)zNT(3103)
对多维函数,可定义与式(370)相同的偏导数。根据该定义可得以下对复矢量的常用导数。aHzz=a*(3104)zHAzz=ATz*(3105)
接下来考虑复参数估计的CRB,假设观测数据为Z=[z0,z1,…,zN-1],观测数据的信号产生模型已知,但其中有L个未知参数,记为θ=[θ1,θ2,…,θL]T。待估计的目标参数为θ的函数,记为p=[p1(θ),p2(θ),…,pK(θ)]T(3106)
记估计子为r(Z),考虑无偏估计子,有∫Ωr(Z)f(Z|θ)dZ=p(3107)
其中Ω表示观测空间,上式两边对θ求偏导可得∫Ωr(Z)f(Z|θ)θHdZ=pθH(3108)上式左边重写为∫Ωr(Z)f(Z|θ)θHdZ=∫Ωr(Z)lnf(Z|θ)θHf(Z|θ)dZ
=Er(Z)lnf(Z|θ)θH(3109)于是pθH=Er(Z)lnf(Z|θ)θH(3110)
简记r(Z)为r,lnf(Z|θ)为lnf,构造如下矢量r
lnfθHH(K L)×1(3111)可算出其自协方差矩阵,并由自协方差矩阵的半正定特性可得Cov(r,r)pθH
pθHHElnfθHHlnfθH≥0(3112)
其中已经用到式(3110)的结论。根据自协方差矩阵特性,下式成立D-pθHI-1(θ)Cov(r,r)pθH
pθHHI(θ)D-pθHI-1(θ)H≥0(3113)
其中D为单位矩阵,进一步将费希尔信息量矩阵定义为I(θ)=ElnfθHHlnfθH(3114)
经简单计算,式(3113)矩阵可化简为Cov(r,r)-pθHI-1(θ)pθHH≥0(3115)
可以得到估计子的误差界为Cov(r,r)≥pθHI-1(θ)pθHH(3116)上式和式(392)实参量估计的CRB是统一的。下面以一个例子介绍CRB的具体计算过程。
例36通信中常用以下观测模型: xn=Asn wn,n=0,1,…,N-1(3117)其中A为复值信道响应,sn为正交调制发送信号,wn是零均值复高斯白噪声序列,方差为σ2,假设发送符号有归一化功率,即sn2=1,求信道估计的CRB。解根据式(295)给出的复高斯概率分布,关于A的一维条件概率分布为p(x|A)=1πσ2exp-xn-Asn2σ2(3118)多变量的联合概率密度为p(x|A)=1πσ2N∏N-1n=0exp-xn-Asn2σ2(3119)基于式(3114)费希尔信息量矩阵的定义及式(350)所示的规范化条件,可以得到与式(339)类似的结论。I(θ)=Elnp(x|θ)θlnp(x|θ)θ=-E2lnp(x|θ)θθ(3120)对数联合概率密度的一、二阶偏导数为lnp(x|A)A=∑N-1n=0snxn-Asnσ2(3121)2lnp(x|A)AA=-1σ2∑N-1n=0sn2=Nσ2(3122)可得到估计信道响应的性能界为CRBA=σ2N(3123)
3.4最大似然估计当被估计量为未知常量时,可以采用比较简单的最大似然估计。最大似然估计可以简便地完成对复杂估计问题的求解,而且,当观测数据足够多时,其性能也是非常好的。因此,最大似然估计在实际中得到了广泛采用。本节只考虑实参数的估计。3.4.1最大似然估计的基本原理设观测矢量x=[x0x1…xN-1]T,待估计参量为θ。观测数据与待估计量的概率依赖关系以概率密度p(x|θ)描述。概率密度是以θ为参量的函数,在观测给定的条件下,如x=x0,p(x=x0|θ)反映了θ取各个值的可能性大小,称p(x=x0|θ)为θ的似然函数,p(x=x0|θ)简记为p(x|θ)。估计问题本质上就是根据观测求出未知量,也就是说,在得到某个观测值x=x0后,如何根据这个观测来确定未知量θ的值。从另外一个角度,似然函数p(x=x0|θ)反映了对不同的θ值观测到x=x0的概率。因为当前的观测结果刚好就是x0,于是可认为实际中观测到x=x0的概率应该尽可能大,即θ取值应该使得p(x=x0|θ)取得最大。使似然函数最大所对应的参数θ作为对θ的估计,称为最大似然估计,记为θ^ml,即θ^ml=argmaxθ p(xθ)或θ^ml=argmaxθ ln p(x|θ)(3124)
下面以两个例子说明最大似然准则的应用。例37设实随机变量满足以下分布: x∝N(0,σ2),有N次独立观测x0,…,xN-1,求σ2的最大似然估计。解单变量似然函数为p(x|σ2)=12πσ21/2exp-x22σ2(3125)多变量联合似然函数为p(x|σ2)=12πσ2N/2exp-12σ2∑N-1n=0x2n(3126)
对数联合似然函数lnp(x|σ2)=-N2ln(2πσ2)-12σ2∑N-1n=0x2n(3127)
令lnp(x|σ2)σ2=0(3128)得到估计结果为σ^2ml=1N∑N-1n=0x2n(3129)
很容易验证2lnp(x|σ2)(σ2)2σ2=1N∑N-1n=0x2n<0(3130)
所以式(3129)求得的是最大值。例38假定观测序列为xn=Acos(2πf0n ) wn,n=0,1,…,N-1(3131)
幅度A和频率f0是已知的,wn是零均值高斯白噪声序列,方差为σ2,求相位的最大似然估计。解似然函数为p(x|)=1(2πσ2)N/2exp-12σ2∑N-1n=0[xn-Acos(2πf0n )]2(3132)
对数似然函数及其导数为lnp(x|)=-N2ln(2πσ2)-12σ2∑N-1n=0[xn-Acos(2πf0n )]2(3133)lnp(x|)=-1σ2∑N-1n=0[xn-Acos(2πf0n )]Asin(2πf0n )(3134)
令上式等于零可得∑N-1n=0xnsin(2πf0n ^ml)=A∑N-1n=0cos(2πf0n ^ml)sin(2πf0n ^ml)
=2A∑N-1n=0sin(4πf0n 2^ml)(3135)
当f0不在0或1/2附近时,上式右边是周期序列求和,当N足够大时,可认为其值近似为零。因此,最大似然估计近似满足∑N-1n=0xnsin(2πf0n ^ml)=0(3136)
展开上式,得∑N-1n=0xnsin2πf0ncos^ml=-∑N-1n=0xncos2πf0nsin^ml(3137)^ml=-arctan∑N-1n=0xnsin2πf0n∑N-1n=0xncos2πf0n(3138)
前面讨论的是单参量的最大似然估计问题,在实际中经常需要同时估计多个参量,例如,在雷达信号的处理中,当检测到目标后,需要同时估计目标的位置、速度等参量,这样的问题称为多参量同时估计,或称为矢量估计。式(3124)所描述的最大似然估计很容易推广到矢量估计的情形。考虑实参数估计问题,设有p个参量θ0,θ1,…,θp-1需要同时估计,定义一个p维矢量θ=[θ0,θ1,…,θp-1]T,那么最大似然估计定义为θ^ml=argmaxθp(x|θ)或者θ^ml=argmaxθlnp(x|θ)(3139)
类似地,最大似然估计的必要条件是p(x|θ)θ|θ=θ^ml=0或者lnp(x|θ)θ|θ=θ^ml=0(3140)
例39设有N次独立观测xn=A wn,其中wn∝N(0,σ2),A,σ2均未知,求A,σ2的最大似然估计。解令θ=[Aσ2]T,则p(x|θ)=12πσ2N/2exp-12σ2∑N-1n=0(xn-A)2(3141)lnp(x|θ)=-N2ln(2πσ2)-12σ2∑N-1n=0(xn-A)2(3142)lnp(x|θ)θ=Nσ21N∑N-1n=0xn-A-N2σ4σ2-1N∑N-1n=0(xn-A)2(3143)
令lnp(x|θ)θθ=θ^ml=0,可求得最大似然估计为θ^ml=A^mlσ^2ml=x-1N∑N-1n=0(xn-x-)2(3144)其中x-表示样本的直接平均,定义为x-=1N∑N-1n=0xn(3145)以下例子给出更一般化线性估计模型中的最大似然估计算法。例310假定待估计矢量为θ=[θ0,θ1,…,θp-1]T,观测数据x=[x0,x1,…,xN-1]T,由如下的线性模型表示x=Hθ w(3146)其中,H是N×p的矩阵,且N>p,H的秩为p,w是多维高斯噪声,其概率模型为w∝N(0,C),其中C为w的协方差矩阵。求θ的最大似然估计。解由式(291)得似然函数为p(x|θ)=1(2π)N/2[det(C)]1/2exp-12(x-Hθ)TC-1(x-Hθ)(3147)要使似然函数最大,只需使下面的表达式最小J(θ)=(x-Hθ)TC-1(x-Hθ)(3148)
基于式(373),式(374),并利用协方差矩阵对称性,将式对θ求导J(θ)θ=-2HTC-1(x-Hθ) (3149)令导数等于零,即可求得最大似然估计,即HTC-1(x-Hθ^ml)=0(3150)解上面的方程式得最大似然估计为θ^ml=(HTC-1H)-1HTC-1x(3151)由上式可以看出,线性观测模型的最大似然估计是观测结果的线性函数,计算简单,因此在噪声协方差矩阵已知情况下是常用的估计算法。3.4.2变换参数的最大似然估计在许多情况下,希望得到估计参量θ的一个函数,例如α=T(θ),如果利用最大似然准则求得θ^,如何进一步求解α^?下面通过几个例子来说明变换参数的最大似然估计的求法。例311设有N次独立观测xn=A wn,n=0,1,…,N-1。其中wn∝N(0,σ2),A为未知参数,σ2已知,求α=eA的最大似然估计。解α的似然函数及其导数为p(x|α)=12πσ2N/2exp-12σ2∑N-1n=0(xn-lnα)2(3152)lnp(x|α)α=1ασ2∑N-1n=0(xn-lnα)=Nασ21N∑N-1n=0xn-lnα(3153)
令上式等于零,可解得α^ml=exp1N∑N-1i=0xi=exp(x-)(3154)由于x-刚好是A的最大似然估计,所以α^ml=exp(A^ml)(3155)可见,α的最大似然估计只需要用A的最大似然估计代入变换式α=eA中就可以求得。因此,如果变换α=T(θ)是一一对应的,那么变换参数后的最大似然估计可以直接由下式得到。α^ml=T(θ^ml)(3156)这一特性称为最大似然估计的不变性。如果变换α=T(θ)不是一一对应的,则不能简单地应用式(3156),下面通过一个例子来加以说明。例312在前例中,假定要求的估计为α=A2,求α的最大似然估计。解由于A=±α,所以本例的变换不是一一对应的,似然函数p(x|α)需要两个概率密度来描述p1(x|α)=12πσ2N/2exp-12σ2∑N-1n=0(xn-α)2,A≥0(3157)p2(x|α)=12πσ2N/2exp-12σ2∑N-1n=0(xn α)2,A<0(3158)
那么,最大似然估计为α^ml=argmaxα {p1(x|α),p2(x|α)}(3159)上式的计算可以分为两步进行。(1) 给定一个α值,如α=α0,比较p1(x|α0)和p2(x|α0)的大小,如果p1(x|α0)≥p2(x|α0)(3160)
则pT(x|α0)=p1(x|α0),对所有的α取值重复以上过程,得到pT(x|α),称pT(x|α)为修正的似然函数。(2) 求出使pT(x|α)最大的α值作为α的最大似然估计,即α^ml=argmaxαpT(x|α)(3161)
对于本例,要比较p1(x|α0)和p2(x|α0)的大小,只需要比较∑N-1n=0(xn-α)2和∑N-1n=0(xn α)2的大小,令J(x,α)=∑N-1n=0(xn-α)2-∑N-1n=0(xn α)2(3162)
上式经整理得J(x,α)=-4∑N-1n=0xnα=-4Nαx-(3163)可见,如果x-≥0,则J(x,α)<0,p1(x|α0)≥p2(x|α0); 如果x-<0,则J(x,α)>0,p2(x|α0)>p1(x|α0),所以修正的似然函数为pT(x|α)=p1(x|α),x-≥0
p2(x|α),x-<0(3164)
当x-≥0时,由pT(x|α)=p1(x|α)可求得α^ml=x-2; 当x-<0时,同样由pT(x|α)=p2(x|α)可求得α^ml=x-2。综合两种情况可得α^ml=x-2(3165)
需要注意的是,如果将A的原始估计A^ml=x-直接代入变换式α=A2中也可以得到式(3165),但这只是一种巧合。参数变换如果不是一一对应的,那么变换后参数的最大似然估计不能简单地通过代入变换式来得到。
3.5贝叶斯估计最大似然估计没有假定被估计量的先验信息,把待估计量看成一个常量。实际环境中待估计量可能是一个随机变量,并且概率密度已知。这时最大似然方法依然可用,但未能充分利用待估计量的信息。这一节介绍贝叶斯估计,是一类能充分利用待估计量概率信息的估计方法。本节只考虑实参数的估计。3.5.1代价函数在估计某个参量θ时,噪声的影响会使估计产生误差,估计误差是要付出代价的,这种代价可以用代价函数来加以描述,记为c(θ,θ^)。一般而言,估计误差Δθ=θ-θ^越小,代价越小。代价函数可表示为c(θ,θ^),c(θ-θ^)或c(Δθ),典型的代价函数有以下几种。(1) 平方代价函数c(θ,θ^)=(θ-θ^)2(3166)
(2) 绝对值代价函数c(θ,θ^)=|θ-θ^|(3167)
(3) 均匀代价函数c(θ,θ^)=1,|θ-θ^|≥ε/2
0,|θ-θ^|(4) 二次型代价函数如果被估计量是p维随机矢量,令误差矢量为Δθ=θ-θ^,则二次型代价函数定义为c(θ,θ^)=‖Δθ‖2S=(θ-θ^)TS(θ-θ^)(3169)式中,‖Δθ‖为误差矢量的范数,S为p×p维的对称非负定的加权矩阵。二次型代价函数实际上是估计误差的加权平方和,当S为单位矩阵时,二次型代价函数是估计误差的平方和。代价函数也可以考虑其他函数,但实际上,最佳估计对代价函数的选择并不敏感。不失一般性,以下考虑单变量,代价函数确定后,可以计算平均代价
c-=E{c(θ,θ^(x))}=∫∞-∞∫∞-∞c(θ,θ^(x))p(θ,x)dθdx(3170)
上式中用θ^(x)代替θ^,强调θ^是从x计算得到,存在函数依赖关系。所谓贝叶斯估计,就是使平均代价最小的估计。式(3170)可以写成c-=∫∞-∞∫∞-∞c(θ,θ^(x))p(θ|x)dθp(x)dx=∫∞-∞c-(θ|x)p(x)dx(3171)
实际应用中p(x)可能是未知的或有复杂的表达式。在一次观测中,如果观测结果为x,则可认为其对应较大的概率。同时由于p(x)是非负的,所以,使平均代价最小近似等价于使条件平均代价c-(θ|x)最小,也就是使c-(θ|x)=∫∞-∞c(θ,θ^(x))p(θ|x)dθ(3172)最小。贝叶斯估计与代价函数的选取有关,使用不同的代价函数将得到不同的贝叶斯估计。3.5.2最小均方误差估计当代价函数为平方代价函数时,平均代价刚好等于估计误差的均方值。先考虑单个参数的估计问题,有
c-=∫∞-∞∫∞-∞[θ-θ^(x)]2p(θ,x)dθdx=E{[θ-θ^(x)]2}(3173)
使平均代价最小等价于使均方误差最小,这时的贝叶斯估计称为最小均方误差估计,记为θ^ms。把平方代价函数代入式(3172),得c-(θ|x)=∫∞-∞[θ-θ^(x)]2p(θ|x)dθ(3174)上式对θ^求导,得c-(θ|x)θ^=-2∫∞-∞[θ-θ^(x)]p(θ|x)dθ(3175)并令导数在θ^=θ^ms处为零,得θ^ms=∫∞-∞θp(θ|x)dθ=E(θ|x)(3176)由于条件平均代价对θ^的二阶导数为正,因此θ^ms所对应的条件平均代价为极小值,即θ^ms为最小均方估计,由式(3176)可以看出,最小均方误差估计为被估计量θ的条件均值。如果被估计量为矢量θ=[θ0,θ1,…,θp-1]T,将二次型代价函数代入式(3170),得c-=∫∞-∞∫∞-∞[θ-θ^(x)]TS[θ-θ^(x)]p(θ,x)dθdx
=E{[θ-θ^(x)]TS[θ-θ^(x)]}(3177)
将二次型代价函数代入式(3172),得c-(θ|x)=∫∞-∞[θ-θ^(x)]TS[θ-θ^(x)]p(θ|x)dθ(3178)
对于θ^(x)求导并令导数在θ^=θ^ms处为零,得-∫∞-∞2S[θ-θ^(x)]p(θ|x)dθ=0(3179)
求解可得θ^ms=∫∞-∞θp(θ|x)dθ=E(θ|x)(3180)
由于2c-(θ|x)θ^θ^T=2S,所以式(3180)求得的是极小值,θ^ms称为最小均方误差估计。3.5.3条件中位数估计当代价函数为绝对值代价函数时,条件平均代价为(考虑单变量情况)c-(θ|x)=∫∞-∞θ-θ^(x)p(θ|x)dθ
=∫θ^(x)-∞[θ^(x)-θ]p(θ|x)dθ ∫∞θ^(x)[θ-θ^(x)]p(θ|x)dθ(3181)
上式对θ^求导,并令导数在θ^=θ^abs处为零,得∫θ^abs-∞p(θ|x)dθ=∫∞θ^absp(θ|x)dθ(3182)由式(3182)可以看出,采用绝对值代价函数的贝叶斯估计θ^abs刚好是条件概率密度p(θ|x)的中位数(Median),所以也称为条件中位数估计,记为θ^med。3.5.4最大后验概率估计当采用均匀代价函数时,条件平均代价为c-=c(θ|x)=∫θ^-Δ/2-∞p(θ|x)dθ ∫∞θ^ Δ/2p(θ|x)dθ
=1-∫θ^ Δ/2θ^-Δ/2p(θ|x)dθ(3183)
很显然,当式(3183)中后面的积分最大时,平均代价最小。要使所述积分达到最大,当p(θ|x)具有类似于高斯分布的特征时,θ^可选为p(θ|x)的最大值所对应的θ(如图32(a)所示),这时对应的贝叶斯估计就是最大后验概率(Maximum a posteriori,MAP)估计,记为θ^map。但是,当p(θ|x)不具备上述特征时,MAP不能得到最优的估计,如图32(b)所示,p(θ|x)具有指数分布的特征,这时θ^map显然不是合理的选择。
图32θ条件概率的不同形式
高斯信号是实际应用中最常见的信号,再加上MAP估计思路直观,实现也比较简便,MAP估计成为一种重要的参数估计准则。下面进一步讨论最大后验概率估计与最大似然估计的联系,最大后验概率估计的数学描述为θ^map=argmaxθ p(θ|x)或θ^map=argmaxθ ln p(θ|x)(3184)
由于p(θ|x)=p(x|θ)p(θ)p(x)(3185)当难以得到p(x)的表达式时,对当前样本,认为p(x)有固定值,这时后验概率密度最大和使p(x|θ)p(θ)最大近似相同,因此,最大后验概率估计为θ^map=argmaxθ p(x|θ)p(θ)或θ^map=argmaxθ[lnp(x|θ) lnp(θ)](3186)上式直观给出了最大后验概率估计和最大似然估计的联系,可以看到最大后验概率估计实际上是在最大似然估计的基础上利用待估计参数的先验知识进行修正。当后验概率密度是可导的函数时,最大后验概率估计的必要条件是p(θ|x)θθ=θ^map=0或lnp(θ|x)θθ=θ^map=0(3187)
式(3187)称为最大后验概率方程。式(3187)所描述的最大后验概率估计很容易推广到矢量估计的情形。假定待估计量为p维矢量θ=[θ0,θ1,…,θp-1]T,后验概率密度为p(θ|x),那么,矢量最大后验概率估计为θ^map=argmaxθ p(θ|x)或θ^map=argmaxθ ln p(θ|x)(3188)
3.5.5贝叶斯估计举例例313设观测为x=A w,其中被估计量A在[-A0,A0]上均匀分布,测量噪声w∝N(0,σ2),求A的最大后验概率估计和最小均方误差估计。解先求最大后验概率估计,因为只有一个观测值,得到条件概率为p(x|A)=12πσexp-(x-A)22σ2
p(A)=12A0,-A0≤A≤A0
0,A>A0(3189)
根据贝叶斯公式p(A|x)=p(x|A)p(A)p(x)(3190)
这时难以得到p(x)的表达式,用下式进行最大后验概率估计Amap=argmaxA[p(x|A)p(A)](3191)
当-A0≤x≤A0时,p(x|A)p(A)的最大值出现在A=x处,所以,A^map=x; 当x>A0时,p(x|A)p(A)的最大值出现在A=A0处,A^map=A0; 当xx,-A0≤x≤A0
A0,x>A0(3192)
易证最大似然估计的结果为A^ml=x。在A满足[-A0,A0]均匀分布的假设下,显然最大后验概率估计的结果更合理。图33给出了不同噪声强度下MAP和ML估计的均方误差性能比较,可见在高噪声强度下,MAP估计性能明显优于ML估计。在实际应用中,如果能获得被估计参数的先验概率信息,可改善估计结果。
图33不同噪声强度下MAP和ML估计性能比较
下面求最小均方估计,由式(3176),得
A^ms=∫∞-∞Ap(A|x)dA=∫∞-∞Ap(x|A)p(A)p(x)dA=∫∞-∞Ap(x|A)p(A)dA∫∞-∞p(x|A)p(A)dA=∫A0-A0A2πσexp-(x-A)22σ2·12A0dA∫A0-A012πσexp-(x-A)22σ2·12A0dA=∫x A0x-A0(x-u)·exp-u22σ2du∫x A0x-A0exp-u22σ2du=x-2σ2∫(z a)/2(z-a)/2u·exp(-u2)duσ∫z az-aexp(-u2/2)du=x-σ{exp[-(z-a)2/2]-exp[-(z a)2/2]}2π[Q(z-a)-Q(z a)](3193)
式中,a=A0/σ,z=x/σ表示归一化观测值,Q(·)为标准正态概率密度函数的概率右尾函数,如下式所示。Q(v)=12π∫ ∞vexp(-x2)dx(3194)例314设有N次独立观测xn=A wn, n=0,1,…,N-1,其中wn∝N(0,σ2),A∝N(μA,σ2A),求A的后验概率估计和最小均方误差估计。解先求后验概率密度p(A|x)=p(x|A)p(A)∫∞-∞p(x|A)p(A)dA=1(2πσ2)N/2exp-12σ2∑N-1n=0(xn-A)212πσ2Aexp-12σ2A(A-μA)2∫∞-∞1(2πσ2)N/2exp-12σ2∑N-1n=0(xn-A)212πσ2Aexp-12σ2A(A-μA)2dA
=exp-121σ2(NA2-2NAx-) 1σ2A(A-μA)2∫∞-∞exp-121σ2(NA2-2NAx-) 1σ2A(A-μA)2dA=exp-12W(A)∫∞-∞exp-12W(A)dA(3195)
式中,x-=1N∑N-1n=0xn为样本均值,且W(A)=1σ2(NA2-2NAx-) 1σ2A(A-μA)2(3196)
注意,式(3195)的分母与A无关,W(A)是A的二次型,经过配方可以把式(3196)写成W(A)=1σ2A|x(A-μA|x)2-μ2A|xσ2A|x μ2Aσ2A(3197)式中,σ2A|x=Nσ2 1σ2A-1(3198)μA|x=Nσ2x- μAσ2Aσ2A|x(3199)
将式(3197)代入式(3195),得p(A|x)=exp-12σ2A|x(A-μA|x)2exp-12μ2Aσ2A-μ2A|xσ2A|x∫∞-∞exp-12σ2A|xA-μA|x2exp-12μ2Aσ2A-μ2A|xσ2A|xdA=12πσ2A|xexp-12σ2A|x(A-μA|x)2(3200)
由式(3200)可以看出,后验概率密度是高斯的。由于最小均方误差估计为被估计量的条件均值,所以A^ms=μA|x=Nσ2x- μAσ2Aσ2A|x=Nσ2x- μAσ2ANσ2 1σ2A(3201)另外,由于最大后验概率估计是使后验概率最大对应的A值,因此,由式(3183)可得A^map=μA|x=A^ms(3202)
3.6线性最小均方误差估计对于随机参量的估计,在3.5.2节中介绍了最小均方误差估计,最小均方误差估计是被估计量的条件均值,这个条件均值通常都是观测的非线性函数,估计器实现起来比较复杂。条件均值的计算需要用到被估计量θ的概率密度p(θ),如果并不知道概率密度p(θ),而只知道θ的一、二阶矩特性,并且希望估计器能用线性系统实现,这时可以采用线性最小均方误差估计。3.6.1随机参量的线性最小均方误差估计线性最小均方误差(Linear Minimum Mean Square Error, LMMSE)估计是一种使均方误差最小的线性估计。假定观测为{xn,n=0,1,…,N-1},那么估计结果定义为观测的线性组合θ^=∑N-1n=0anxn b(3203)估计的均方误差为Mse(θ^)=E[(θ-θ^)2]=Eθ-∑N-1n=0anxn-b2(3204)
选择一组最佳系数an和b,使式(3204)的均方误差达到最小。注意,引入系数b是因为观测x和被估计量θ的均值可能不为零。如果观测x和被估计量θ的均值都为零,那么系数b就可以省略。均方误差对系数b求导,并令导数等于零,得Mse(θ^)b=-2Eθ-∑N-1n=0anxn-b=0(3205)
经整理后得b=E(θ)-∑N-1n=0anE(xn)(3206)将上式代入式(3204),得Mse(θ^)=E∑N-1n=0an[xn-E(xn)]-[θ-E(θ)]2(3207)令a=[a0,a1,…,aN-1]T,x=[x0,x1,…,xN-1]T,那么,式(3203)可表示为θ^=aTx b(3208)
而式(3206)可改写为b=E(θ)-aTE(x)(3209)估计的均方误差可表示为Mse(θ^)=E{[aT(x-E(x))-(θ-E(θ))]2}
=E[aT(x-E(x))(x-E(x))Ta]-E[aT(x-E(x))(θ-E(θ))]-
E[(θ-E(θ))(x-E(x))Ta] E[(θ-E(θ))2]
=aTCxa-aTcθx-cTθxa Cθ(3210)
其中,Cx=E{[x-E(x)][x-E(x)]T}为观测的协方差矩阵,cθx=E{[θ-E(θ)][x-E(x)]}是待估计量θ与观测x的协方差矢量,且Cθ是θ的方差。均方误差对a求导,得Mse(θ^)a=2Cxa-2cθx(3211)
令导数等于零可解得a=C-1xcθx(3212)
将式(3212)和式(3209)代入式(3208),得θ^lmmse=cTθxC-1xx E(θ)-cTθxC-1xE(x)(3213)
上式经整理后得到线性最小均方估计为θ^lmmse=E(θ) cTθxC-1x[x-E(x)](3214)如果θ和x的均值为零,则θ^lmmse=cTθxC-1xx(3215)将式(3215)代入式(3210)中可以得到最小均方误差的表达式Mse(θ^lmmse)=Cθ-cTθxC-1xcθx(3216)线性最小均方误差能得到简洁的估计子,并且可定量分析估计误差,因此应用广泛。例315设有N次独立观测数据xn=A wn,n=0,1,…,N-1(3217)其中,A在(-A0,A0)上服从均匀分布,wn是零均值高斯白噪声,方差为σ2,且A与wn统计独立,求A的线性最小均方误差估计。解可以把观测模型写成矢量形式x=AeN×1 w(3218)其中,eN×1是N维的全1矢量。根据题意,E(A)=0,σ2A=E(A2)=(2A0)2/12,E(x)=0且Cx=E(xxT)=E[(Ae w)(Ae w)T]=E(A2)eeT σ2I=σ2AeeT σ2I(3219)
cAx=E(Ax)=E[A(Ae w)]=E(A2)e=σ2Ae(3220)
其中I为单位矩阵,由式(3215),得A^lmmse=cTAxC-1xx=σ2AeT[eeTσ2A σ2I]-1x=σ2Aσ2eT[eeTσ2Aσ2 I]-1x(3221)
根据矩阵求逆中常用的ShermanMorrison公式(B cdT)-1=B-1-B-1cdTB-11 dTB-1c(3222)式(3221)经整理后,得A^lmmse=σ2Aσ2A σ2/Nx-=A20/3A20/3 σ2/Nx-(3223)
3.6.2线性最小均方误差估计的几何解释可以用向量空间的概念来解释线性最小均方误差估计。在以下的讨论中假定被估计量θ和观测x都是零均值的,如果不是零均值,则总可以定义零均值的随机变量θ′=θ-E(θ),x′=x-E(x)。
考虑一个如图34所示的由所有随机变量构成的线性空间,θ和x0,x1,…,xN-1可以看作线性空间的点(每个点由一个向量描述,如图34(a)所示)。定义两个向量的内积为(x,y)=E(xy),向量的长度为随机变量方差的均方根,即‖x‖=E(x)2,如果(x,y)=E(xy)=0,称这两个向量是正交的。由于随机变量是零均值的,因此,两个向量正交和随机变量不相关等价。线性最小均方误差估计的目的是用x0,x1,…,xN-1的线性组合来逼近θ。先考虑只有一个样本的情况,即用θ^=a0x0来逼近θ,易见最优的估计就是θ到x0的投影(如图34(a)所示),这时误差为θ到x0的垂直距离,从几何的角度易知这时误差向量的长度最短。当使用两个样本x0,x1时,θ^=a0x0 a1x1,注意到x0,x1的线性组合张成一个平面,这时最优估计就是θ到该平面的投影(如图34(b)所示),同时误差向量ε正交于由{x0,x1}所张成的平面,易证这时误差向量的长度最短。
图34线性最小均方误差估计的几何解释
一般情况下,线性估计子用x0,x1,…,xN-1的线性组合来近似θ,相当于在x0,x1,…,xN-1所张成的子空间上找一个点,使其与θ的距离最小。根据误差定义Mse(θ^)=E[(θ-θ^)2]=Eθ-∑N-1i=0aixi2=‖θ-∑N-1i=0aixi‖2(3224)
可见最小均方误差等价于误差向量ε=θ-θ^长度的平方最小。根据上面的论述,这时误差向量将正交于x0,x1,…,xN-1所张成的子空间,即ε⊥span{x0,x1,…,xN-1}(3225)根据正交的定义,可以得出E[(θ-θ^)xn]=0,n=0,1,…,N-1(3226)上式是线性最小均方估计非常重要的正交原理,也可由式(3224)直接对加权系数求解得到。由正交原理可知,线性最小均方误差估计可以通过使估计误差和每个观测数据正交得到。利用正交原理,加权系数很容易求出Eθ-∑N-1n=0anxnxj=0,j=0,1,…,N-1(3227)或者写成∑N-1n=0anE(xnxj)=E(θxj),j=0,1,…,N-1(3228)用矩阵表示为E(x20)E(x0x1)…E(x0xN-1)E(x1x0)E(x21)…E(x1xN-1)E(xN-1×0)E(xN-1×1)…E(x2N-1)a0a1aN-1=E(θx0)E(θx1)E(θxN-1)(3229)
或者写成Cxa=cθx(3230)
所以a=C-1xcθx(3231)θ的线性最小均方估计为θ^lmmse=aTx=cTθxC-1xx(3232)上式与式(3215)一致。下面考虑随机矢量的线性最小均方误差估计。设待估计量为θ=[θ0,θ1,…,θp-1]T,观测矢量x=[x0,x1,…,xN-1]T,线性估计可表示为θ^=Ax b(3233)
式中,b是p×1维的矢量,A是p×N维的矩阵,所有估计的均方误差和可表示为Mse(θ^)=E∑p-1i=0(θi-θ^i)2=E[(θ-θ^)T(θ-θ^)](3234)将式(3233)代入式(3234),得Mse(θ^)=E[(θ-Ax-b)T(θ-Ax-b)](3235)分别对A和b求导并令导数等于零,得Mse(θ^)b=-2E[θ-Ax-b]=0(3236)Mse(θ^)A=-2E{[θ-Ax-b]xT}=0(3237)
上式中对矩阵的求导与式(372)相同,定义为对矩阵每个元素求偏导数,结果为矩阵,且大小与A相同。由上述两式,可得b=E(θ)-AE(x)(3238)A=CθxC-1x(3239)式中,Cθx=Cov(θ,x)=E{[θ-E(θ)][x-E(x)]T}(3240)Cx=Cov(x)=E{[x-E(x)][x-E(x)]T}(3241)
由此可得线性最小均方误差估计为θ^lmmse=E(θ) C-1xCθx[x-E(x)](3242)式(3242)在表达形式上除θ是矢量外与式(3214)相同。线性最小均方误差估计的性质总结如下(这里不做具体证明): (1) 线性最小均方误差估计是无偏估计。(2) 线性最小均方误差估计的均方误差阵要小于任何其他线性估计的均方误差阵,即E[(θ-θ^lmmse)(θ-θ^lmmse)T]≤E[(θ-θ^)(θ-θ^)T](3243)式中,θ^是任意的线性估计。上式中符号≤作用于矩阵时,含义如下: 如果A≥B,则A-B是半正定矩阵。式(3243)证明从略。(3) 对于线性变换α=Hθ b,线性最小均方误差估计具有不变性,即α^lmmse=Hθ^lmmse b(3244)
(4) 线性最小均方误差估计具有叠加性,即如果α=θ1 θ2,则α^lmmse=θ^1,lmmse θ^2,lmmse(3245)3.7最小二乘估计在前面介绍的几种估计方法中,最小均方误差估计、最大后验概率估计需要知道被估计量的先验概率密度,最大似然估计需要知道似然函数,线性最小均方误差估计需要知道被估计量的一、二阶矩,如果这些概率密度或矩未知,就不能采用这些方法,这时可以采用最小二乘(Least Square)估计。最小二乘估计对统计特性没有做任何假定,因此,它的应用非常广泛。在线性最小均方误差估计中,选择均方误差作为估计量性能好坏标准的度量,也就意味着估计与真实参数之间的差别平均来说达到最小。而在最小二乘估计中,则是使给定的观测数据与信号或者无噪声的观测数据之差的平方和达到最小。采用例310的线性观测模型,待估计量为θ=[θ0,θ1,…,θp-1]T,观测信号用矢量和矩阵可表示为x=Hθ w(3246)
式中,x=[x0,x1,…,xN-1]T,w=[w0,w1,…,wN-1]TH=h00h01…h0(p-1)h10h11…h1(p-1)h(N-1)0h(N-1)1…h(N-1)(p-1)
估计偏差的平方和可表示为J(θ^)=[x-Hθ^]T[x-Hθ^]=∑N-1n=0xn-∑p-1j=0hnjθ^j2(3247)
在没有任何先验信息情况下,可以认为最优估计就是使得上述的平方和误差最小。当待估计参数的重要性存在差别时,可以对每个参数的误差进行加权,观测与估计偏差的加权平方和可表示为Jw(θ^)=[x-Hθ^]TW[x-Hθ^]=∑N-1j=0∑N-1n=0xn-∑p-1k=0hnkθ^kWnjxn-∑p-1k=0hjkθ^k(3248)其中W为加权矩阵。最小二乘估计就是使J(θ^)最小的估计,记为θ^ls; 加权最小二乘(Weighted Least Square)估计则是使Jw(θ^)最小的估计,记为θ^lsw。接下来考虑如何求解,将式(3247)重写为J(θ^)=xTx-xTHθ^-θ^THTx θ^THTHθ^(3249)根据式(373)和式(374)的矢量函数求导公式,得到xTxθ^=0(3250)xTHθ^θ^=HTx(3251)θ^THTxθ^=HTx(3252)θ^THTHθ^θ^=2HTHθ^(3253)
求J(θ^)对θ^的导数,并令导数等于零,得J(θ^)θ^=-2HTx 2HTHθ^=0(3254)由此可解得最小二乘估计为θ^ls=(HTH)-1HTx(3255)求Jw(θ^)对θ^的导数,并令导数等于零,得Jw(θ^)θ^=-2HTW[x-Hθ^]=0(3256)由此可解得加权最小二乘估计为θ^lsw=(HTWH)-1HTWx(3257)最小二乘估计具有如下特点: (1) 对于线性的观测模型,最小二乘估计和加权最小二乘估计都是线性估计,对测量噪声的统计没有做任何假定,应用十分广泛。(2) 当测量噪声的均值为零时,即E(wi)=0时,最小二乘估计和加权最小二乘估计都是无偏估计。(3) 对于加权最小二乘估计,如果有一些模型的知识,如E(w)=0,E[wwT]=R,当W=R-1时,估计误差的方差达到最小,这时加权最小二乘估计和式(3151)的最大似然估计是一致的。3.8信号检测基础信号检测的目的是在干扰、噪声背景下判断目标信号是否存在。当存在多种类型目标信号时,还需要判断信号类型。根据目标信号模型的不同分确定性信号检测和随机信号检测进行讨论,本节只讨论实信号。3.8.1确定性信号检测这里以确定性信号为例介绍信号检测的基本概念和方法。假设目标信号sn的波形已知,得到了接收信号xn的N个样本,在此基础上判决其中是否包含目标信号,可定义如下两种事件(或称假设)。H0:xn=wn,n=0,1,…,N-1
H1:xn=sn wn,n=0,1,…,N-1(3258)
接下来考虑如何对上述事件进行检验,其过程可以看成是对观测空间划分成不同判决区域,观测样本x=[x0,x1,…,xN-1]落在那个区域就判决成相应事件。一般可基于事件的概率来划分判决区域,定义两种事件的后验概率分别为p(H1|x),p(H0|x),后验概率给出了在得到观测样本后,对H0、H1事件发生概率的认知。一种直观的判决准则是,在观测样本的基础上,哪个事件发生的概率大就认为发生了该事件,即当p(H1|x)>p(H0|x)时,认为H1事件(假设)为真。反之若p(H1|x)
p(H0|x)的观测空间。具体的数学描述如下式所示。H1:p(H1|x)p(H0|x)>1(3259)
上式中的后验概率一般难以直接计算。基于贝叶斯公式有p(H1|x)p(x)=p(x|H1)p(H1)(3260)式(3259)等价于H1:L(x)=p(x|H1)p(x|H0)≥γ(3261)
其中p(x|H0),p(x|H1)是3.4节最大似然估计中提到的似然函数,L(x)称为似然比(Likelihood Ratio,LR),有时为了计算便利,使用其对数形式ln[L(x)],称为对数似然比(Logarithm Likelihood Ratio,LLR)。γ=p(H0)/p(H1)表示判决门限。
计算似然函数需要知道xn的概率密度,实际中可通过假设噪声符合某种概率分布实现,这里以最常见的假设: wn是零均值、方差为σ2的高斯白噪声为例,分析如何具体应用式(3261)进行检验。根据假设可得到如下似然函数: p(x|H1)=1(2πσ2)N2exp-12σ2∑N-1n=0(xn-sn)2(3262)p(x|H0)=1(2πσ2)N2exp-12σ2∑N-1n=0x2n(3263)
利用上式可计算如下似然比和对数似然比: L(x)=exp-12σ2∑N-1n=0(xn-sn)2-x2n(3264)l(x)=lnL(x)=-12σ2∑N-1n=0(xn-sn)2-x2n(3265)采用与式(3261)等价的检验准则: lnL(x)>lnγ,经简单计算,可得式(3261)等价于对下式中T(x)的检验。H1:T(x)=∑N-1n=0xnsn=xTs>γ-(3266)其中γ-=σ2lnγ 12∑N-1n=0s2n(3267)注意到T(x)可以看成是匹配滤波器的输出,确定性信号的匹配滤波器定义如下: h(n)=s(N-1-n),n=0,1,…,N-1
0,其他(3268)
匹配滤波器的输出为y(n)=x(n)h(n),并且可验证T(x)=y(N-1)。在已知噪声概率特性情况下,上述分析均可用,下面将上述推导拓展到高斯有色噪声下的信号检测问题。记x=[x0,x1,…,xN-1]T,s=[s0,s1,…,sN-1]T,w=[w0,w1,…,wN-1]T,检测问题可描述为H0:x=w,w~N(0,C)
H1:x=s w(3269)
判决准则重写为H1:L(x)=p(x|H1)p(x|H0)>γ(3270)可得到似然函数,似然函数比及对数似然比为p(x|H1)=1(2π)N2det12(C)exp-12(x-s)TC-1(x-s)(3271)p(x|H0)=1(2π)N2det12(C)exp-12xTC-1x(3272)l(x)=lnL(x)=-12(x-s)TC-1(x-s)-xTC-1x=xTC-1s-12sTC-1s(3273)
利用lnL(x)>lnγ准则进行判决,得到H1:T(x)=xTC-1s>γ 12sTC-1s=γ-(3274)同理,可以由C-1s构造滤波器,而T(x)是滤波器某个时刻的输出,这种滤波器称为广义匹配滤波器。T(x)>γ-定义了观测空间的一种划分规则,对所有满足T(x)>γ-的观测x计算积分可以得到H1成功检测的概率。根据事件的不同,可定义4种不同概率,其定义和计算方式如表31所示。
表31二元信号检测中的事件和概率计算
事件定义标号H1为真,判H1成立H1检测概率P(H1|H1)=∫T(x)>γ-p(x|H1)dx
H0为真,判H0成立H0检测概率1-P(H1|H0)
续表
事件定义标号H0为真,判H1成立虚警概率P(H1|H0)=∫T(x)>γ-p(x|H0)dx
H1为真,判H0成立漏检概率1-P(H1|H1)
在实际应用中上表不同的事件会有不同的代价,例如在雷达信号处理中,漏检造成的代价远大于虚警概率。考虑更一般化的情况,为表中四种不同事件定义各自的代价因子,并定义如下平均代价: C=[C00P(H0|H0) C01P(H1|H0)]P(H0)
[C10P(H0|H1) C11P(H1|H1)]P(H1) (3275)其中C00,C01,C10,C11分别为各事件的代价因子,经简单计算可得
C=[C00(1-P(H1|H0)) C01P(H1|H0)]P(H0)
[C10(1-P(H1|H1)) C11P(H1|H1)]P(H1)
=C00P(H0) (C01-C00)P(H0)P(H1|H0)
C10P(H1)-(C10-C11)P(H1)P(H1|H1)
=C00P(H0) C10P(H1)
∫x∈Z1(C01-C00)P(H0)P(x|H0)-(C10-C11)P(H1)P(x|H1)dx(3276)
其中Z1表示判决H1成立的判决区域。可见要使代价最小,在判决区域内的积分函数取值为负,即在判决区域内满足(C01-C00)P(H0)P(x|H0)(C01-C00)P(H0)(C10-C11)P(H1)(3278)
上式称为信号检测的贝叶斯准则。从上式可以看到,如果认为正确的检测不产生额外代价,即式(3275)中C00=C11=0,同时虚警和漏检有相同的代价,即式(3275)中C01=C10=1,则式(3278)的检测准则等价于式(3261)。
以上的贝叶斯准则需要确定代价因子和先验概率,实际中往往很困难,特别是先验概率的确定。实际中往往对虚警概率和漏检概率有明确的要求,可以从虚警概率和漏检概率的约束考虑信号检测,一种常见的准则是在虚警概率恒定的约束下最小化漏检概率,称为纽曼-皮儿森准则(NeymanPearson,NP)。令PF,PM分别表示虚警概率和漏检概率,α表示允许的虚警概率上限。Z0,Z1分别表示信号检测的判决区域,问题模型可表示为minZ0,Z1PMPF=α(3279)
即在满足虚警概率要求条件下,寻找观测空间的划分准则,使得漏检概率尽量小。可利用拉格朗日乘子法求解上述问题,构造如下目标函数J=PM λ(PF-α)(3280)其中λ为拉格朗日乘子。目的是要确定一种观测空间的划分,使J最小,将虚警概率和漏检概率的计算表达式代入,得到J=∫Z0p(x|H1)dx λ∫Z1p(x|H0)dx-α
=λ(1-α) ∫Z0[p(x|H1)-λp(x|H0)]dx(3281)
上式中第一项非负,要使J尽可能小,则上式积分结果尽可能小,一种直观的做法是将所有使得p(x|H1)-λp(x|H0)<0的观测空间归入Z0,或将使得p(x|H1)>λp(x|H0)成立的观测空间定义为Z1,从而得到如下H1的检测准则: H1:p(x|H1)p(x|H0)>λ(3282)
上式中λ为未知量,不能直接应用。λ可根据约束条件确定,但难以得到闭式解。一般可采用迭代方式求解,即从λ的一个初始值开始,由式(3282)得到观测空间的分割Z0,Z1。再根据如下约束条件进行验证,如果满足则停止迭代,否则调整λ取值重新验证。∫Z1p(x|H0)dx=PF(3283)
3.8.2随机信号检测随机信号检测所采用的准则与确定性信号一致。但由于对随机信号只能知道其统计特征,因此具体计算有所不同,下面看一个简单例子,考虑以下两种假设。H0:x=w,w~N(0,σ2I)
H1:x=s w(3284)
其中sn是平稳随机信号,先考虑一个简单情况,sn是高斯信号,均值为零,方差为σ2s,以上两种假设可描述为H0:x~N(0,σ2I)
H1:x~N(0,(σ2s σ2)I)(3285)采用似然比及对数似然比进行判断,得到如下H1的判决准则: H1:L(x)=p(x|H1)p(x|H0)>γ(3286)
得到似然函数及似然比如下: P(x|H1)=1(2π(σ2s σ2))N2exp-12(σ2s σ2)∑N-1n=0x2n(3287)P(x|H0)=1(2πσ2)N2exp-12σ2∑N-1n=0x2n(3288)l(x)=lnL(x)=N2lnσ2σ2s σ2 12σ2sσ2(σ2s σ2)∑N-1n=0x2n(3289)
采用l(x)>lnγ判决准则,对上式略加整理,可得到如下等效判决准则: H1:T(x)=∑N-1n=0x2n=xTx>lnγ-N2lnσ2σ2s σ22σ2(σ2s σ2)σ2s=γ-(3290)
可以看到H1的判决依赖于采样信号的能量,当能量大于门限值时,判断存在目标信号,因此上述检测方法又称为能量检测器。能量检测器可推广到更一般的情况,即检测信号具有任意的协方差矩阵Cs,考虑以下两种假设H0:x~N(0,σ2I)
H1:x~N(0,Cs σ2I)(3291)得到似然函数及似然比如下: P(x|H1)=1(2π)N2det12(Cs σ2I)exp-12xT(Cs σ2I)-1x(3292)P(x|H0)=1(2πσ2)N2exp-12σ2xTx(3293)l(x)=lnL(x)=12ln(σ2)Ndet(Cs σ2I)-12xT[(Cs σ2I)-1-1σ2I]x(3294)
可得到如下等价于l(x)>lnγ的判决准则: H1:T(x)=σ2xT1σ2I-(Cs σ2I)-1x>γ-(3295)利用矩阵求逆定理(A BCD)-1=A-1-A-1B(DA-1B C-1)-1DA-1(3296)上式中令A=σ2I,B=D=I,C=Cs,有
(σ2I Cs)-1=1σ2I-1σ41σ2I C-1s-1(3297)
代入式(3295),得到如下判决准则:
H1: T(x)=xT1σ2C-1s 1σ2I-1x=xT(σ2C-1s I)-1x
=xT[(σ2I Cs)C-1s]-1x=xTCs(σ2I Cs)-1x>γ-(3298)
和式(3290)相比,上式相当于在能量检测器基础上引入加权矩阵W=Cs(σ2I Cs)-1,得到加权检测量T(x)=xTWx,称为广义能量检测器。注意到式(3298)可写成H1:T(x)=xTs^>γ-
s^=Cs(σ2I Cs)-1x(3299)
信号检测的过程可以看成是先得到目标信号的一个估计,然后再利用估计得到的信号构造一个如式(3266)所示的匹配滤波器。注意到上式中的s^实际是目标信号的线性最小均方误差(LMMSE)估计。说明如下: 根据3.6节的结论,基于观测结果x,目标参数θ的线性MMSE估计为θ^=CθxC-1xx(3300)将上式用于从x中估计s可得Csx=E(sxT)=E(ssT swT)=E(ssT)=Cs(3301)Cx=E(xxT)=E(ssT swT wsT wwT)
=E(ssT) E(wwT)
=Cs σ2I(3302)
上式中假设目标信号和噪声不相关,这在实际中一般可以满足。以上讨论都基于假设在H0和H1条件下的似然函数是可知的,但实际情况经常没这么理想。一个最简单的例子是描述噪声特性的参数,如方差,经常是未知的。目标信号的参数也有可能是未知的,例如在无源雷达和声呐系统中,接收机对来波信号的频率是未知的。估计存在未知参数时,一般可考虑两种方法。一种是先解决H0和H1条件下的未知参数估计问题,例如用最大似然准则θ^0=maxθ0p(x|θ0,H0)θ^1=maxθ1p(x|θ1,H1)(3303)
然后将估计结果用于似然比准则进行判决H1:L(x)=p(x|θ^1,H1)p(x|θ^0,H0)>γ(3304)
这种方法称为广义似然比准则。如果先验概率分布已知,则可以使用贝叶斯准则计算似然比H1:L(x)=p(x|H1)p(x|H0)=∫p(x|θ1,H1)dθ1∫p(x|θ0,H0)dθ0>γ(3305)
下面以一个例子进一步说明,基于式(3290)且假设目标信号具有非零均值H0:x~N(0,σ2I)
H1:x~N(ue,Cs σ2I)(3306)其中e是N×1全1矢量,得到似然函数及似然比如下: P(x|H1)=1(2π)N2det12(Cs σ2I)exp-12(x-ue)T(Cs σ2I)-1(x-ue)(3307)P(x|H0)=1(2πσ2)N2exp-12σ2xTx(3308)
u的最大似然估计为u^=1NxT(Cs σ2I)-1e(3309)
和式(3294)的推导相似,可得到以下判决准则T(x)=-(x-u^e)T(Cs σ2I)-1(x-u^e) 1σ2xTx>γ-
xT[(Cs σ2I)-1-1σ2I]x Nu^2>γ-
xTCs(σ2I Cs)-1x Nu^2>γ-(3310)
和式(3294)相比,上式的不同之处在于引入了非零均值的影响。上述讨论可推广到多种信号的检测,考虑如下模型: H0:x=w,w~N(0,C)
Hi:x=si w,i=1,…,K(3311)其中K表示有多少种不同的信号。为了实现最低的检测错误概率,可采用最大后验概率检测准则Hi:p(Hi|x)>p(Hj|x),j=0,1,…,K,j≠i(3312)假设各种假设具有相等概率,可采用最大似然检测准则Hi:p(x|Hi)>p(x|Hj),j=0,1,…,K,j≠i(3313)
在实际应用中经常先计算各种假设对H0的似然比Li(x)=p(x|Hi)p(x|H0),i=1,2,…,K(3314)
如果所有Li(x)的取值均小于1,则判决H0为真,否则判决具有最大似然比的假设为真。特别要指出的是,上述方法只能应用于各种假设具有相等概率的情况,如果概率不相等,则必须直接基于式(3312)计算。本章习题1. 假定观测x0,x1,…,xN-1是独立同分布的,且在(0,θ)上服从均匀分布,求θ的最大似然估计。2. 证明式(339)I(θ)=Elnp(x|θ)θ2=-E2lnp(x|θ)θ2
3. 设N次观测为xi=A wi,i=0,1,…,N-1,其中A为未知的确定信号,噪声wi (i=0,1,…,N-1) 相互独立并服从相同的分布N(0,σ2),噪声的方差已知。(1) 求A的最大似然估计。(2) 是否为无偏估计?(3) 是否为有效估计?估计的方差等于多少?4. 在习题3中,如果A和σ2都是未知的,试分别求其最大似然估计,并讨论它们的无偏性。5. 设观测量为xi=A wi,其中wi~N(0,σ2),如果A为未知常量,求A的最大后验概率估计和最小均方误差估计。6. 设有N次独立同分布的观测x0,x1,…,xN-1,且单个样本在(0,θ)内服从均匀分布,证明:
Elnp(x|θ)θ≠0,θ>07. 平坦信道的信道估计问题可以用下式描述xi=Asi wi,i=0,1,…,N-1其中,wi是零均值高斯白噪声,方差为σ2,si是已知的训练序列,xi是接收信号,求估计A的CRB。证明有效估计量存在,并求它的方差。当N→∞时估计的方差会怎样变化?8. 设有两次观测x0=As0 w0x1=As1 w1其中,w=[w0,w1]T是零均值高斯随机矢量,协方差矩阵为C=σ21ρρ1
求估计A的CRB,并将它与w为高斯白噪声的情况(即ρ=0)进行比较,解释当ρ=±1时会怎样?9. 考虑高斯噪声下的单频频率估计问题xi=Aejωn φ wi,i=0,1,…,N-1其中wi是复高斯白噪声,方差为σ2。A,ω,分别为幅度,频率和相位,求(1) A、为已知参数时频率估计的CRB; (2) A、ω、都是未知参数时频率估计的CRB。10. 设观测的概率密度为p(x|θ),θ为未知常量,待估计量为α=g(θ),证明: 估计α的CRB为Var(α)≥g(θ)θ2-E2lnp(x|θ)θ2
11. 假定α=g(θ)=Aθ b,假设θ-是θ的有效估计,即Var(θ-)=I-1(θ),证明: α=Aθ b也是有效估计量,即对于线性变换,估计量的有效性得以保持。12. 考虑一个正弦参数的估计问题,假定观测为xi=Acos(2πf0i φ) wi,i=0,1,…,N-1
其中,wi为零均值高斯白噪声,方差为σ2,A、f0、为未知参量,且A>0,00
0,A<0其中λ>0,wi是零均值,方差为σ2n的高斯白噪声,与A相互独立,试求A的最大后验概率估计。17. 对无偏估计子,最小均方误差等价于最小方差。对线性估计子θ^=∑N-1i=0aixi
引入约束条件E[xi]=siθ,则最小方差无偏估计子参数如下: a=C-1ssHC-1s
其中C是xi的协方差矩阵,证明所得估计子是无偏估计子并求其估计误差。18. 设x(t)=AcosΩ0t w(t),通过采样对幅度A作线性估计。在Ω0t=0,Ω0t=π/3处采样两次,记为x0,x1。并设E[A]=0,E[A2]=2,E[w0]=E[w1]=0,E[w20]=E[w21]=1,E[w0w1]=0,E[Aw0]=E[Aw1]=0,求: (1) 线性最小均方误差估计A^lms; (2) 线性最小均方误差估计的误差。19. 对于x=Hθ w的实系数线性观测模型,假定测量噪声的均值为零,证明最小二乘估计误差的均方误差阵为Pθls=E(θ-θls)(θ-θls)T=(HTH)HTRH(HTH)-1
其中R=E[wwT]。20. 观测信号模型为yi=α βxi wi,i=0,1,…,N-1。通过采样对幅度a作线性估计,其中xi已知,wi为未知干扰,求α,β的最小二乘估计。
21. 设有零均值实平稳随机过程x(t),t是[0,T]内的一点。若已知x(0)和x(T),利用x(0)和x(T)求x(t)的最佳线性估计。22. 设有两种如下假设,观测次数为N次H0,zk=nk
H1,zk=2 nk(k=1,2,…,N)
其中nk是满足零均值,方差为σ2n的正态分布,假设p(H0)=p(H1)=0.5,求: (1) 最小错误概率准则下的判决表达式; (2) 虚警概率和检测概率(结果由误差函数表示)。
23. 设有两种如下假设,观测次数为N次H0,zk=nk
H1,zk=1 nk(k=1,2,…,N)
其中nk是满足零均值,方差为σ2n的正态分布,试构造一个虚警概率为0.1的NP检测器,求相应的检测概率。24. 考虑如下四元数字通信系统: H0,zk=nk
H1,zk=1 nk
H2,zk=2 nk
H3,zk=3 nk
其中nk是满足零均值,方差为σ2n的正态分布,各个假设的先验概率相等。试设计一个四元信号的最佳检测系统,由每一个接收符号zk检测其发送符号。
评论
还没有评论。