误差估计

2024-08-28

误差估计（共8篇）

误差估计篇1

1 前言

电力系统状态估计是电力系统调度、控制、安全评估等应用的基础, 也是能量管理系统 (EMS) 的核心组成部分[1]。状态估计的根本目的是通过估计计算排除量测数据中的噪声干扰, 获得对系统当前状态的最优估计值, 因此, 量测数据中的噪声特性对估计的效果有着显著的影响。目前广泛采用的状态估计大多是基于基本加权最小二乘法及其改进方法[2,3,4,5]。一般近似认为量测数据中只含有零均值白噪声或高斯白噪声, 即量测误差服从标准正态分布, 基于这种假设的加权最小二乘估计是无偏估计。然而, 实际中的量测误差并不完全符合上述假设, 误差中除了含有服从正态分布的随机误差外, 往往还包括系统误差。系统误差是指误差的数值是固定的或按非统计规律变化的误差[6]。系统误差通常是由测量传感器未校准或温度、湿度、电磁环境变化使传感器测量精度下降造成的。与随机误差相比, 系统误差的数值要大的多, 并且在一定时间内持续存在, 因此, 对量测精度的影响更大。

目前, 对于系统误差的排除一般是通过硬件校准方法实现的, 即对传感器进行校准, 然而, 在校准前并不知道那些传感器需要校准, 因此只能采取逐个排查的方法, 费时、费力并且效率较低。与硬件校核方法相对, 可以采用一种软件校核方法, 即通过对传感器采集数据特性的分析, 识别出存在系统误差的量测, 并且通过特定的算法从量测数据中剔除系统误差, 实现量测数据的自动校准。

2 系统误差对状态估计结果的影响

在给定网络接线、支路参数和量测数据的条件下, 电力系统的量测方程为:

z=h (x) +v (1)

其中, z为量测向量, x为系统状态变量向量, h (x) 为x的非线性函数, v为误差向量。理论上, 人们通常假设v中只含有随机误差, 并且服从均值为0的标准正态分布, 即E (v) =0。实际上, v除了包含随机误差外, 有时还含有系统误差, 即v=se+ne, se为系统误差, ne为随机误差。根据误差的性质有E (se) =uk≠0, E (ne) =0, 因此有

给定测量向量z以后, 状态估计向量undefined是使目标函数 (3)

达到最小值的x。h (x) 通常为x的非线性函数, 可以通过迭代法求解。当求出undefined后, 可以用h (undefined) 代替h (x) , 最终的残差为

其中,

对 (4) 两端取数学期望, 即

由 (1) 式, 并用h (undefined) 代替h (x) , 则, undefined, 即undefined, 由于H (undefined) , 所以undefined。因此当量测误差中含有系统误差时, 按照加权最小二乘原则得到的估计值是真实状态的有偏估计值。通过上述分析可知, 系统误差的存在将会使估计值偏离真实值, 从而影响状态估计的精度。

3 系统误差的检测

由于系统误差是由传感器未校准造成的, 所以与随机误差不同, 系统误差一般是局部的, 即只在某些量测中存在系统误差。对系统误差进行检测, 识别出含有系统误差的量测是处理系统误差的第一步。可用一种简单实用的残差检测法用于识别系统误差, 即在多次状态估计计算之后通过对残差的检测来判断量测中是否含有系统误差。如果量测误差服从均值为0的标准正态分布, 则估计得到的残差, 如果出现即可认为量测数据中含有系统误差。此残差检验法充分利用已有状态估计的计算结果, 不需要增加新的计算模块, 只要记录若干次残差并进行简单计算即可判断出含有系统误差的量测。在实际应用中可以将作为判断系统误差的依据, 其中M为门槛值, 也就是说, 如果某一量测残差的数学期望超过门槛值即可认为该量测中含有系统误差。

4 计及系统误差的状态估计方法

4.1 量测校准模型

发现系统误差的根本目的是消除系统误差, 目前还没有一种统一的、通用的、完善的系统误差处理方法, 如果想直接计算出系统误差的大小是比较困难的, 但是通过分析系统误差的特点不难发现, 系统误差可以分为线性和非线性两大类, 由于这两种误差的存在, 使得系统状态真值与测量值之间近似存在某种线性或非线性的函数关系, 如果能够辨识出这种函数关系就可以求出剔除系统误差的量测值, 达到消除系统误差影响的目的。基于上述思想提出的状态估计的数学模型如下:

假设剔除系统误差后的量测值 (校准量测) 与含有系统误差的量测值 (未校准量测) 之间存在如下函数关系:

其中, undefined为包含系统误差的量测向量 (未校准量测) , z为经校准后, 不含系统误差的量测向量 (校准量测) p为待辨识的参数, f为校准函数

将 (1) 式和 (7) 式展开为一阶泰勒级数得

其中,

将 (8) 代入 (9) 得

undefined

其中, v′=Fz·v, 将 (10) 式写成矩阵形式

其中, ,

通过求解式 (11) 中的未知向量Δy可以得到待辨识的参数p, 进而得到校准函数f, 通过求解校准函数的反函数即可得到剔除系统误差的量测值 (12)

Zc为剔除系统误差后的量测向量

以上便是量测校准的完整过程, 求解向量是整个过程的关键, 式 (11) 给出的是基于一个断面数据的某一量测的计算公式, 通常情况下一个时间断面数据的冗余度不足以计算多个量测的, 可以采用一定时间窗口内的多个断面的数据增加冗余度。扩展后的状态向量为:

其中为选定的时间窗内的数据断面个数。相应的, 式 (11) 扩展为:

利用最小二乘法求解 (15) 式, 可得参数p。

4.2 具有量测校准功能的状态估计方法

通过上面的分析可知, 当量测数据中含有系统误差时, 按照传统加权最小二乘法得到的估计值将偏离真实值, 通过参数辨识的方法可以得到量测校准函数, 利用校准函数剔除系统误差的影响, 达到提高状态估计精度的目的, 具体的应用流程如下:

1) 记录M次常规状态估计的残差向量Vi=[v1 v2 …vn], i=1, 2, …M, n为量测个数, 计算量测残差的均值 (数学期望) 向量

将均值大于门槛值的残差对应的量测作为待校准量测;

2) 选择校准函数f, 对于线性系统误差可以令undefined, 对于非线性系统误差可以令undefined;

3) 记录K个时间断面的量测数据, 这里需要注意的是K的取值不能太大, 以免计算量过大;

4) 初始化状态向量 (xundefined, i=1, 2, ...K) 和校准参数向量 (p0) , 并使迭代变量j=1;

5) 计算Fzij, Fundefined, Hundefined, i=1...K以及每次迭代过程中的雅克比矩阵J (Hj) , 其中表示某一时间断面, j为记录迭代次数的变量;

6) 计算每次迭代的状态估计值undefined, 并计算相应的量测误差

(16)

7) 求解最小二乘状态估计问题, 计算状态变量和参数p的修正值 (Δxundefined, i=1...K) 和 (Δpj) , 并且更新状态变量和参数p的值xundefined=xundefined+Δxundefined, pj=pj-1+Δpj;

8) 判断是否收敛, 即|Δxundefined|<ε且|Δpj|<ε, 如果收敛转至97, 否则更新迭代变量j=j+1并转至57;

9) 得到参数p并代入 (12) 计算出剔除系统误差后的量测值, 完成一轮校准。

10) 返回到1) 并用校准后的量测代替原来的未校准量测, 计算M次常规状态估计的残差, 检测是否存在未校准量测, 若存在则重复2) —9) , 直到所有存在系统误差的量测都被校准。

5 结论

以上提出的计及系统误差的状态估计算法, 是与传统的状态估计方法相结合的, 可以方便地加以应用, 提高状态估计精度。

参考文献

[1]于尔铿.电力系统状态估计[M].北京:水利电力出版社, 1985.

[2]张海波, 李林川.电力系统状态估计的混合不良数据检测方法[J].电网技术, 2001, 25 (10) :17-20.ZHANG Hai-bo, LI Lin-chuan.A hybrid ap-proach for detection of bad data in power system stateestimation[J].Power System Technology, 2001, 25 (10) :17-20.

[3]周良松, 赵卫东.一种改进的Hachtel状态估计方法[J].电力系统及其自动化学报, 2006, 18 (5) :42～45.ZHOU Liang-song, ZHAO Wei-dong.Modi-fied Hachtel Method for state estimation[J].Pro-ceedings of the EPSA, 2006, 18 (5) :42～45.

[4]钱峰, 龚庆武.基于IGG法的电力系统状态估计[J].电力系统自动化, 2005, 29 (3) :36～41.QIAN Feng, GONG Qing-wu.State estimationof power system based on IGG method[J].Auto-mation of Electric Power System, 2005, 29 (3) :36～41.

[5]李强, 周京阳, 于尔铿等.基于混合量测的电力系统状态估计混合算法[J].电力系统自动化, 2004, 29 (19) :31～36.LI Qiang, ZHOU Jing-yang, YU Er-keng, LIU Shu-chun, WANG Lie.A Hybrid Algorithm forPower System State Estimation Based on PMU Meas-urement and SCADA Measurement[J].Automa-tion of Electric Power System, 2004, 29 (19) :31～36.

[6]梁晋文.误差理论与数据处理[M].中国计量出版社, 2001.

误差估计篇2

天线方向性误差对DOA估计影响及其校正

天线阵元的方向性误差会导致基于特征值分解的DOA估计算法性能的`下降,因此有必要对阵元方向性误差进行校正.首先建立了存在方向性误差情况下DOA估计的数学分析模型,提出利用空间平滑法和改进空间平滑法对方向性误差进行校正.仿真实验表明,方向性误差对DOA估计的影响主要表现在谱峰峰值的降低.空间平滑法和改进空间平滑法均可以有效地对方向性误差进行校正.

作者：于斌尹成友黄冶罗怀中 YU Bin YING Cheng-you HUANG Ye LU Huai-zhong 作者单位：解放军电子工程学院,合肥,230037 刊名：电光与控制 ISTIC PKU英文刊名：ELECTRONICS OPTICS & CONTROL 年，卷(期)： 14(3) 分类号：V271.4 TN91 关键词：误差 DOA 空间平滑校正

误差估计篇3

解析函数边值问题在奇异积分方程方面有广泛的应用, 它们在弹性力学、流体力学方面也有重要的应用.在实际问题中边界曲线发生摄动的情况经常遇到的, 比如空间的弹性基本问题和平面的、空间的断裂力学问题.所以, 不少学者在致力于边值理论分析的同时, 也考虑了将理论运用到解决实际问题中。因此, 研究解析函数边值问题关于边界曲线的稳定性有很大的实际意义。

1 问题的提出

之差, 本文的主要结果将给出这两个奇异积分之差的估计式。

2 预备知识

则由引理1可得

3 主要结果

证明:.考虑

则文献[3]定理4可得

一方面, 由推论1可得

另一方面, 由Xω+, ∏, g的Hölder连续性可得

这里, υ∈ (0, 1) 是任意取定的数, 利用上述两方面的估计式可得

则

综上所述

定理的结论得证。

摘要：为了得到带根号的Riemann边值问题边值问题关于边界曲线摄动的的稳定性, 因此本文讨论了与之相对应的一类奇异积分关于积分曲线摄动的误差估计。

关键词：带根号Riemann边值问题,奇异积分,摄动,误差估计

参考文献

[1]路见可.解析函数边值问题[M].武汉:武汉大学出版社, 2009:413-415.

[2]曾乔, 林峰.一类奇异积分关于积分曲线摄动的误差估计[J].四川师范大学学报 (自然科学版) , 2015 (01) .

误差估计篇4

多用户MIMO技术( Multi-user MIMO,MU-MIMO) 可以获得多址容量的增益[1]。其实现方法有:脏纸编码( Dirty Paper Code,DPC) 被证明能够达到理论的信道容量[2],但实现复杂度极高; 迫零波束成型( Zero-Forcing Beam Forming,ZFBF)[3,4]实现简单,在用户数充分大时,多用户分集效应使ZFBF也可以获得接近DPC的总速率性能,因此得到了广泛的研究。

经典ZFBF算法假设基站可以获得理想的信道的状态信息( Channel State Information,CSI) ,而在实际的系统中,只能将信道信息量化后反馈给基站。由于频谱资源有限,在信道相干时间内反馈信道仅仅可以提供有限的反馈比特数,因此会带来速率损失。文献[5 - 8]对有限反馈的量化误差对总速率性能带来的影响进行了分析,并提出了改进方案。但都没有关注信道估计误差对系统总速率性能的影响。

本文对多用户MIMO下行链路的ZFBF算法进行了分析,将信道估计误差考虑在内,分有用户调度与无用户调度2 种场景分析信道估计误差对ZFBF总速率性能的影响,并在此基础上以最大化有效传输速率为目标,计算得到最优的训练序列长度。

1 系统模型

考虑如图1 所示的单小区多用户MIMO通信系统的下行链路,假设小区中有K个用户,基站配有M根天线,每个用户均配备单天线。用户k的信道用1 × M的行向量hk表示,则系统的信号模型为:

式中,H = [h1H,…,hKH]H为K × M维的信道矩阵;x为M × 1 维的发射信号向量; n为K × 1 维的高斯白噪声向量,均值为0,方差为1; H和n的元素相互独立且服从复高斯分布: CN(0,1) ; y为K × 1 维的接收信号向量。

对第k个用户,记sk、wk和Pk分别为对应的发送符号向量、预编码波束向量和发射功率缩放系数,则用户k对应的发射信号为,接收信号为:

式中,第1 项是期望信号; 第2 项为其他用户信号造成的干扰; 第3 项为高斯噪声。ZFBF算法的基本思想是通过预编码向量的设计使期望接收端的干扰信号为0,即选择wj,使hkwj= 0,j ≠ k 。

2 信道估计误差对总速率性能的影响

在多用户MIMO系统中,下行波束的构造很大程度上依赖于信道的状态信息,信道状态信息的获取需要由发送端发送训练序列,并在接收端进行信道估计。假设训练序列长度为T ,用M × T维的矩阵S表示训练序列矩阵,接收端采用最小均方误差( Minimum Mean Square Error,MMSE) 算法进行估计,由于在MIMO信道中采用酉阵作为训练矩阵是MSE最优的[9],因此按照下面的标准选取S :

式中,P为基站的发射功率,对应的接收信号为:

信道矩阵的MMSE估计为[9,10]:

根据式( 5) ,信道矩阵可以表示为估计矩阵与误差矩阵的和:,估计矩阵和误差矩阵E相互独立并且它们的元素均为零均值、循环对称的复高斯分布,方差分别为:

ZFBF实现的前提条件是同时传输的用户数不大于发射天线数,即K ≤ M 。当K > M时需要先进行用户调度,即从所有用户中选出M个用户,然后在选出的用户组内构造波束矢量。下面分2 种情况进行分析。

① K ≤ M

当K ≤ M时,不需要进行用户调度。基站根据信道估计值构造迫零单位波束矢量,即

假设对所有用户等功率分配,则每个用户的接收信号为:

式( 8) 括号内的第2 项表示多用户干扰。式中,ek表示用户k的信道估计误差向量。接收机将干扰项视为加性高斯噪声,则对应的遍历总速率为:

记wk为假设不存在信道估计误差时得到的迫零波束向量,则对应的遍历总速率为:

定义为理想情况下的总速率与存在信道估计误差时的总速率差值,即信道估计误差带来是速率损失,可以得出这一差值的上界为[9]:

由这个上界可以得到M个用户的总速率的近似值:

另外,从式( 11) 可以看出,随着P的增长,速率损失逐渐增大,因此,在高信噪比条件下要维持系统的性能需要提高信道估计的精度。所幸的是,当P → ∞时,ΔRk的上界将趋于一个常数,即

② K > M

当用户数很大时,采用穷搜索实现用户调度的复杂度过高,因此经常使用低复杂度的次优用户选择算法。半正交用户选择( Semiorthogonal User Selection,SUS) 算法[3]是在所有用户中选出一组半正交( 空间分离度高) 且信道质量好的用户,可以证明基于SUS算法的ZFBF是渐进最优的,本文中使用这一算法进行用户调度,再构造波束向量。由于SUS算法与信道估计过程相互独立,所以SUS算法对 ΔRk造成的影响相对较小,因此仍然可以用式( 11) 作为 ΔRk的近似值。定义如下总速率的近似表达式:

式中,。当K充分大时,RZF可以由M log2(1 + P log2K / M) 来近似,得到总速率的另一个近似表达式:

3 最优训练长度

假设相干时间长度为Tc,用于信道估计的训练时隙数为T,( T ≥ M) ,每个时隙长度为 τ 。则用于数据传输的时间为Tc- Tτ ( 忽略信息反馈、估计有效信道信息等操作) 。假设实际的数据传输速率为R ,令 α = τ / Tc,可以定义有效传输速率为,定义最优训练时隙数为使有效传输速率最大的T ,即

考虑信道估计误差带来的影响,使用RAPPROX2来近似R 。因此,近似的最优训练时隙数为:

求( 1 - αT) RAPPROX2对T的导数并令其为0,可以求出的值。由于求导后的表达式难以直接求出闭合解,当 α 取值较小且P取值较大时,对式( 17) 中的求优表达式做进一步近似:

接下来对T求导并令其值为0,可得TOPT的另一个近似值:

下面将通过仿真分析和的近似程度,以及TOPT随P的变化规律。

4 仿真结果

为了验证上述分析的准确性,将根据文中数学模型得到的分析结果与实际系统仿真进行比对。仿真建立在图1 所示的系统模型下,发射天线数目M= 4 ,接收天线数目为1,信道为平坦瑞利衰落,其元素服从均值为0、方差为1 的复高斯分布。

T = 6 、无用户选择时,总速率- 平均SNR的仿真图如图2 所示。图2 中“理想CSIT”是不考虑信道估计误差的理想情况下系统总传输速率曲线,“实际速率”是采用MMSE信道估计的实际仿真结果,“近似速率”则根据式( 12) 对传输速率的近似表达得出的仿真结果。可以看出近似速率与实际速率的变化趋势基本相同,且差值并不大,尤其是在高SNR时,近似速率非常接近实际速率,充分说明了近似式( 12) 的准确性。

K = 100,T = 6 时的总速率- 平均SNR的仿真图如图3 所示。图3 中“近似速率1”和“近似速率2”分别对应RAPPROX1和RAPPROX2。从图3 中可以看出,高SNR时RAPPROX1与实际速率非常接近,说明采用SUS算法进行用户调度对 ΔRk的界定没有什么影响。RAPPROX2的近似效果略差,这是由于使用M log2( 1 + P log2K / M) 来近似RZF带来的影响,但是差值并不大,在信噪比30 d B时,RAPPROX2与实际速率的差值大约为1 bps/Hz。因此,在一些定性分析的场景中完全可以使用RAPPROX2来代替实际速率。

K = 100 ,P = 15 d B时,有效总速率- 训练时隙数的仿真图如图4 所示。α 取值分别为0. 005、0.01 和0. 02。图4 中的实线为实际的仿真值,虚线则是根据RAPPROX2得到的近似值。从图4 中可以看出,虽然近似值和实际值在数值上有比较明显的差异,但是它们的变化趋势很相似,极值点的位置也基本相同,因此根据RAPPROX2得到的能较好地刻画TOPT。另外,注意到TOPT附近的几个T值对应的有效速率与最大的有效速率值很接近,如果考虑到仿真中存在的误差,取这几个值作为最优值都是可以接受的,考虑实际系统设计时用于信道估计的开销越小越好,在实际应用中训练序列长度可以在TOPT附近选取较小的值。

K = 100,α 分别取0. 001 和0. 005 时,最优训练时隙数- 平均SNR的仿真图如图5和图6 所示。图5和图6 中的实线为TOPT的实际仿真值,“近似1”和“近似2 ”分别代表和。可以看出和的近似程度都较好。因此利用式( 19)可以近似分析TOPT随 α 、P 、K和M等参数的变化规律。式( 19) 也显示出,随着信噪比的增加,需要的训练序列长度减小,这与图中的变化趋势也是一致的。

5 结束语

本文在MU-MIMO环境下,对存在信道估计误差时的ZFBF传输过程建立模型,分析了信道估计误差带来的速率损失,并得出下行总速率的近似表达; 进而以最大化有效总速率为目标,求得训练序列的最优长度; 最后通过仿真对前面的分析进行了验证。仿真结果表明,基于文中模型的分析结果与实际情况吻合度较高,能够较为准确地刻画信道估计误差对ZFBF传输速率性能的影响,通过这种近似,可以更方便地求得训练序列的最优长度,为系统的设计提供参考。

参考文献

[1]GESBERT D,KOUNTOURIS M,HEATH J R W,et al.From Single User to Multiuser Communications:Shifting the MIMO Paradigm[J].IEEE Sig.Proc.Magazine,2007,24(5):36-46.

[2]COSTA M.Writing on Dirty Paper[J].IEEE Transactions on Information Theory,1983,29(3):439-441.

[3]YOO T,GOLDSMITH A.The Optimality of Multiantenna Broadcast Scheduling Using Zero-forcing Beamforming[J].IEEE Journal on Selected Areas in Communications,2006,24(3):528-541.

[4]赵辰,刘应状,朱光喜.VLST系统中ZF检测算法的研究[J].无线电通信技术,2006,33(2):37-39.

[5]CAIRE Giuseppe,JINDAL Nihar,KOBAYASHI Mari,et al.Multiuser MIMO Achievable Rates with Downlink Training and Channel State Feedback[J].IEEE Transactions on Information Theory,2010,56(6):2 845-2 866.

[6]SHEN Hong,XU Wei,JIN Shi,et al.Joint Transmit and Receive Beamforming for Multiuser MIMO Downlinks With Channel Uncertainty[J].IEEE Transactions on Vehicular Technology,2014,63(5):2 319-2 335.

[7]KIM Haksoo,YU Heejung,LEE Yong H.Limited Feedback for Multi-cell Zero-forcing Coordinated Beamforming in Time-varying Channels[J].IEEE Transactions on Vehicular Technology,2014,64(6):2 349-2 360.

[8]NI Wei,CHEN Zhuo,SUZUKI Hajime,et al.On the Performance of Semi-orthogonal User Selection with Limited Feedback[J].IEEE Communications Letters,2011,15(12):1 359-1 361.

[9]MARZETTAT L.BLAST Training:Estimating Channel Characteristics for High Capacity Space-Time Wireless[C]∥Preceedings of the 37th Annu.Allerton Conf.Communication,Control,Computing,Monticello,1999:22-24.

误差估计篇5

多元传感器融合处理能通过综合多角度、多类传感器的不同观测, 实现数据率、精度和发现概率的提高。要实现这些指标的提升, 需要对多个传感器的观测值进行综合估计。但是由于不同传感器的数据率可能并不相同, 或者即便相同数据率的传感器对同一目标的观测时间也不一样, 这就使得数据融合首先要解决数据的时间配准[1]。

数据融合按照时间配准的方法可分为同步融合和异步融合两类。同步融合[2]是将多个传感器的观测值同步到一个时刻, 并进行综合参数估计;异步融合[3]并不需要做时间同步, 而是直接对每一个观测进行单独的参数估计过程。因此异步融合的数据率通常比同步融合更高, 对不同传感器的观测时序不作任何假设或要求。

1多源数据融合的一般方法

多源传感器的数据融合一般分为三个阶段:配准、参数估计和滤波, 而一般单一传感器的目标跟踪不需要配准和参数估计。因此就工程实践来说, 多源传感器的目标跟踪与单传感器的跟踪在软件实现步骤上基本一致, 只不过多源融合程序增加了时间配准和参数估计的过程。整个的融合过程如图1所示。

1.1融合架构

信息融合方法从架构上可分为分布式和集中式两种[4,5]。分布式融合[6]是对每一个传感器观测进行独立的跟踪, 产生局部航迹, 之后对所有的局部航迹进行融合, 产生全局航迹。这样做的优点是融合中心得到的信息更多 (如速度、加速度等) , 有利于时间配准中对观测数据在时间维度上的插值和外推;缺点是由于需要对单独传感器进行独立跟踪, 因此难以通过增加传感器提高发现概率。

集中式融合将各传感器的原始观测值进行统一融合处理, 直接产生全局航迹, 这样做的好处是能提高系统整体的发现概率, 当单一传感器发现概率达不到系统要求时, 通过数据融合使全局发现概率达到检测要求。

1.2时间配准

为了对多个传感器数据进行融合, 必须将多个传感器的数据和融合历史数据配准到同一个时刻, 然后才能进行参数估计。时间配准包含两个任务:1.选择合适的时刻;2.对数据进行预测或内插。

一般情况下, 选定某一个主传感器的观测时刻为融合时刻, 或者使用融合系统自身处理节拍的时刻作为融合时刻, 所有其他传感器的数据统一地通过预测或内插调整至融合时刻。这种方法称为同步融合, 好处是能够在每一个时刻得到更多的观测数据, 从而进行精确的估计;坏处是融合数据率低于传感器的原始数据率, 并且更多的传感器无法提高融合数据率, 而只能提高精度。

直接对传感器的观测数据进行融合[7], 并在每次收到任何一个传感器观测数据后都将该数据与历史数据做实时融合并给出融合结果。这种方法称为异步融合, 好处是更多的传感器可以得到更高的数据率的同时也能提高融合精度, 坏处是实时观测数据只能和历史数据进行融合, 而非最新的观测数据, 从而影响融合精度。

2基于误差概率的异步融合

本文方法以集中式融合为架构, 将多传感器的观测值直接与全局航迹关联和融合, 实现全局航迹的实时更新, 是一种集中式异步融合方法。

2.1集中式异步融合处理

集中式异步融合处理是将传感器观测与全局航迹的外推值进行融合, 因此每一次融合分为三个步骤:1) 全局航迹外推;2) 参数估计;3) 航迹滤波。

本文使用的航迹滤波方法采用匀加速运动卡尔曼滤波。实际上, 基于误差概率的异步融合方法与采用的滤波算法没有关系, 更好的、更符合目标运动特征的滤波算法能够得到更好的结果。

匀加速运动模型的卡尔曼滤波算法[5]如下:

若要将航迹位置外推至时刻Tk+1, 则预测位置

2.2误差概率估计

传感器通常会上报多个维度的观测数据, 大部分是以传感器为原点的极坐标。不同类型的传感器其距离和方位误差差别巨大, 例如基于电磁波的传统雷达通常距离精度高, 方位精度差;而基于可见光的激光雷达距离精度很差, 但方位精度极高。可以对这些不同传感器分别建立误差模型, 并以此模型为基础进行数据融合[8]。

误差概率估计的融合方法如下:

1) 首先假设传感器数据在每一维度的误差Ei是已经建模好的;

2) 然后将误差转换为对应的概率Pi, 但所有维度的概率必须保证归一化, 即∑Pi≡1;

3) 最后利用这个误差概率将观测值与全局航迹外推值进行融合, 即zi=yiPi+xi (1-Pi) , 其中x代表观测值, y代表全局航迹外推值, z代表融合值, i代表维度。

例如以传统二维雷达为例, 假设其理论距离误差为100米, 方位误差为3度, 最大量程为300公里。由于卡尔曼滤波是在直角坐标系下进行的, 因此距离和方位误差会相互耦合, 那么最大距离误差概率为Rmax×3×π/180, 方位误差仍然约为3度。可由下面的方法计算得到距离和方位的误差概率:

由于跟踪是在直角坐标系下进行, 因此在数据融合过程中, 是将极坐标下的观测值转换到直角坐标系中, 与全局航迹的直角坐标外推值进行融合, 因此还需要将极坐标下的误差概率转换到直角坐标系中。直角坐标系下的误差概率是原始观测方位值的函数, 即:

显然Px和Py也都是归一化的。

需要注意的是在跟踪开始时刻, 由于航迹需要时间收敛至接近真值, 在未收敛之前的这段时间采用全局航迹进行误差概率估计会导致结果不准确。解决方法是在跟踪初始阶段直接采用观测值进行滤波, 相当于将这一阶段的观测值误差概率设为0。

3仿真实验

3.1实验方法

仿真实验模拟了两个传统两坐标雷达之间的数据融合。首先按照三次样条曲线等间隔生成了1000组观测值 (以目标匀速曲线运动为条件) , 然后根据两个雷达的数据率分别在这1000组数据中抽取若干作为各雷达的观测, 并加入方位误差 (距离误差可忽略不计) 。产生的数据如图2所示。

对产生的数据分别采用直接异步融合、同步估计融合、异步估计融合三种方法进行数据融合。其中直接异步融合是将两个雷达的观测点当作单传感器直接进行跟踪;同步估计融合将其中数据率低的雷达时间配准到数据率高的雷达观测时刻, 然后进行同步融合;异步估计融合就是本文上一节讲到的方法。同步估计融合和异步估计融合全部采用误差概率估计方法对观测值进行参数估计。

3.2实验结果

上述三种方法都采用相同的匀加速卡尔曼滤波器, 最终的结果如图3、图4和表1所示。

4总结

通过仿真对比试验, 本文所提出的基于误差估计的集中式异步融合算法即保证的数据率, 又使得精度得到了提高, 其原因是同步融合中对观测值的外推时间跨度更大, 因此引入的误差比采用全局航迹的外推更大, 特别是在目标转弯时表现非常明显。

摘要：为了提高多传感器数据融合中的精度, 本文提出了基于误差估计的异步融合处理方法 , 该方法将传感器观测值与融合航迹的预测值直接融合, 同时保证多传感器融合结果的数据率和融合之后的精度。本文进行了对比仿真实验, 结果表明该方法在直线和转弯条件下的估计精度都比同步方法更高。

关键词：异步融合,航迹跟踪,误差概率

参考文献

[1]刘钊.多传感器时空一致及其信息融合技术研究[D].国防科学技术大学, 2008.

[2]刘钊, 戴斌.多传感器信息融合中时间同步方法的研究[J].计算机仿真, 2009, 26 (6) :124-127.

[3]何俊, 张华, 刘继忠.基于耦合HMM的多传感异步信息融合分类策略[J].计算机应用研究, 2009 (8) :2836-2838.

[4]宋文尧.卡尔曼滤波[M].科学出版社, 1991.

[5]韩崇昭.多源信息融合[M].清华大学出版社, 2006.

[6]陶贵丽, 邓自立.广义系统多传感器分布式融合降阶Kalman滤波器[J].科学技术与工程, 2006, 6 (6) :661-668.

[7]刘磊, 林雪原, 周旻.多传感器组合导航系统的多尺度异步信息融合算法[J].海军航空工程学院学报, 2012 (3) :276-280.

误差估计篇6

关键词：有效值,卷积窗加权算法,不同步误差,量化误差

以微处理器为基础的电工测量仪表在计算电力系统有效值时, 信号频率的不稳定会带来不同步采样误差[1,2]。许多学者就克服这种测量误差进行了研究, 如小波分析法[3], 三点法[4], 微分采样法[5], 各方法都有其优点。文献[6]提出电气参量卷积窗加权算法, 可基本消除不同步采样的测量误差, 且权函数有固定解析式, 算法简单便于实现。但对于加权算法计算有效值时, ADC带来的量化误差问题却没有进行详细的研究。本文分析了不同步采样误差与相对频偏和窗阶数的关系, 重点推导了卷积窗下有效值量化相对误差关于采样频率和量化位数的解析式。

1电气参量卷积窗加权算法分析

1.1有效值定义

设交流电电压信号u (t) 的周期为T, 有效值定义为: 在一个信号周期内, 通过某纯阻负载所产生的热量与一个直流电压在同一负载上产生的热量相等时, 该直流电压的数值就是交流电压的有效值, 通常也称真有效值, 即

对于标准正弦信号, 其有效值为:

在电信号采样中, 即使等间隔采样的采样周期与信号周期存在整数倍关系, 但由于信号频率不稳等原因, 导致用式 (1) 计算有效值时, 使得信号周期T只能用一个近似值T0代入, 会出现不同步误差。

1.2有效值的卷积窗加权算法[6]

采样法计算有效值的实质是对信号加长度为T0, 中心在t0的矩形窗函数:

其中w1为矩形函数:

k阶卷积窗加权有效值测量形式为:

其中为k阶卷积窗。

w1 (t) 为宽度为T0的矩形窗函数, 故wk (t) 的时域总宽度为k T0, 其傅里叶变换为:

1.3加权算法的不同步采样误差分析

对于最高谐波次数为M, 基频为f1的周期信号u2 (t) , 可以表示为:

式中:A0为信号u2 (t) 的直流分量 (即u2 (t) 的平均值) ; Am及 φm分别为第m次谐波的幅值和初相位。

定义加窗信号h (t, t0) 为:

Uk2 (t) 为加窗信号h (t, t0) 的平均值, 根据傅里叶变换的定义, 它显然等于h (t, t0) 的直流分量H (0, t0) , H ( f, t0) 为h (t, t0) 的傅里叶变换, 则:

根据定义可知A0= U2, 式 (8) 第二项便是k阶加权有效值平方Uk2 (t0) 在非同步采样 (T0≠ T1) 时所造成的误差。

为了进一步分析误差的大小, 引入相对同步偏差x , 即T0= (1 + x) T1, f1= (1 + x) f0。

在电力系统中, 由于基频的波动范围很小, 一般来说x ≪ 1且最高谐波次数M < 20, 故mx ≪ 1 , 这时对式 (8) 进行化简:

可以看出测量误差正比于x (1 + x) 的k次方, 比例系数为窗函数起始处信号对其均值的偏离量, 若选择采样信号段起始处使其采样值刚好等于信号的均值, 则测量误差将等于0。

若测量值随t0的不同而起伏变化, 其数学期望和方差分别为:

因此, 由非同步采样造成的测量相对误差 σ 为:

由于x ≪ 1 (0.01量级) , 可以看出随着权阶k的增大, 不同步采样的测量误差急剧减小。对于电力系统中信号, 一般只需k取2或3, 便可把误差降低到次要的、可以忽略的程度。

2量化误差分析

由于利用上述算法基本上可消除不同步采样对各电气参量所造成的测量误差, 因此信号的量化对测量造成的误差将变得重要。在分析有效值卷积窗加权算法的量化误差之前, 先对信号进行时间上的离散[7]。

对以任意时刻t0为中心, 周期为T1的信号x (t) 进行等间隔采样, 假设N为一个估计周期T0中的采样点数, 那么采样间隔 Δt = T0N。采样周期为采样窗阶数k, 采样序列包含k N个数据, 采样序列记为xn, 相应的时域权重为ws:

其中n = 0, 1, ⋯, k N - 1, wk为k阶卷积窗, 由文献[6] 决定。

对于电压有效值U的加权离散表达式为:

设ADC量化位数为p比特, 对电压un量化得到ûn, 其差值就是量化误差, 定义为:

设U为电压信号u (t) 的真有效值, Û为在同步采样条件下由于量化而得到的有效值, 则根据式 (15) 有:

设u (t) 的最大峰峰值为Um, 且为ADC的满度值, 由于采样起点的随机性, 每个采样点的量化误差en可以认为在Um2-p[-1 2, 1 2] 服从均匀分布, 则:

不同采样点的量化误差可以看成是相互独立无关的, 则有效值的期望和方差为:

其中:与信号的波形有微弱的相关性, 但主要依赖于卷积窗的阶数k, 对于正弦信号[8]有: A1= 100, A2= 1.15, A3= 1.29, A4= 1.39, A5= 1.47 。因此由于量化而造成的测量误差eq为:

由于打点计数器的采样频率fs是可控制的, 对于量化相对误差eq与量化位数p及采样频率fs之间的关系, 可表示为:

通过式 (21) 可以看出, 采样数据量化引起的测量相对误差正比于量化位数而反比于采样频率的开方根, 而增大卷积窗阶数k会使Ak也增大, 故只使误差略微减小。

从噪声角度分析[9], 如果ADC的量化误差噪声近似为白噪声, 则有:

式中:L为过采样比p为量化噪声功率;X m为ADC的满幅度值。

可以看出, 采样频率在奈奎斯特频率的基础上每增加一倍, 相当于量化位数增加1 2位, 这与式 (21) 得出的结论相吻合。同时可知, 若原信号不变, 量化信噪比每增加6 d B, 测量误差eq减小一倍。

3数值模拟

对一般正弦电压信号u (t) = sin (2π (1 + x) f0t + φ) 进行Matlab数值仿真。电力系统中信号频率波动范围约为49.5 ∼ 50.5 Hz, 取T0= 20 ms ( f0= 50 Hz) , 则相对频偏| x|< 0.01 。通过式 (12) 可给出当卷积窗阶数k分别为1, 2, 3时, 相对频偏与不同步测量误差的关系, 如图1所示。

图2中实线表示量化位数p为8, 采样频率fs= 3 200, 卷积窗阶数k = 3时, 通过式 (21) 算出的量化相对误差与不同步误差组成的综合误差 (两误差平方和的开方) 与相对频偏之间关系的估计值, 实际上在此条件下量化误差为主要误差。空心点表示相同条件下的实际仿真结果, 其过程为:在一个估计周期内取360组不同的中心点tj, 量化每组un为ûn, 分别代入式 (14) 求出360组不同的有效值Û (tj) , 求得它们的相对误差 (标准差、均值) 。各点误差幅度波动较大是由于不同的相对频偏导致每采样周期量化误差叠加程度不同, 在同步采样 (x = 0) 时, 各周期量化误差完全叠加, 因此误差最大。

表1列出了在正常有效值为220 V时, 用本文理论算出的综合绝对误差大小 (相对频偏为1%) 。可以看出:

(1) 估计值和仿真值基本吻合;

(2) 当k = 1时, 不同步误差为主要误差, 增大采样频率及量化位数对综合误差影响较小;

当k ≥2时, 量化误差为主要误差, 增大k使综合误差略微减小。

(3) 量化位数较低时, 增大采样频率, 可使量化绝对误差由10-1数量级减少为10-2甚至更低的数量级。

4结论

随着卷积窗阶数k的增加, 不同步误差以相对频偏xk速度减少。各级卷积窗解析式[6]已有, 因此算法比较简单。

减少有效值的量化误差主要是增大ADC的量化位数p及采样频率fs, 增大卷积窗阶数k会使误差略微减小。采样频率每增加一倍, 相当于量化器位数增加1 2位。若给定误差eq范围及p值, 通过式 (21) 可确定采样频率fs。

综上, 电力系统中有效值的测量, 通过卷积窗加权算法及适当的采样频率便可把ADC带来的综合误差控制在可忽略水平。

参考文献

[1]李沂乘.采样计算方法测量交流电压有效值误差分析[J].电子测量与仪器学报, 2008 (z1) :62-65.

[2]KAMPIK M, LAIZ H, KLONZ M.Comparison of three accurate methods to measure ac voltage at low frequencies[C]//Proceedings of the 16th IEEE Instrumentation and Measurement Technology Conference.Venice, Italy:IMTC, 1999:429-433.

[3]郑势, 王秋光, 石雨涛.采用小波变换的有效值和频率测量方法[J].电子应用, 2007, 26 (2) :101-103.

[4]葛庆庆, 邵天章.交流电信号有效值高精度测量算法研究[J].电测与仪表学报, 2012, 49 (11) :11-14.

[5]赵明旭, 周佩娟.交流模拟量参数精确测量仪表设计[J].东北电力技术, 2011 (6) :36-37.

[6]张介秋, 陈砚圃, 梁昌洪, 等.交流电电气参量高精度测量的加权算法[J].仪器仪表学报, 2005, 26 (2) :172-176.

[7]WANG Yang, ZHANG Hong-ming, WU Qing-wei, et al.Improvement of photonic ADC based on phase-shifted optical quantization by using additional modulators[J].IEEE Journals&Magazines, 2012, 24 (7) :566-568.

[8]张介秋, 梁昌洪, 陈砚圃.卷积窗及其在电力系统参量估计中的应用[J].电子学报, 2004, 32 (12) :2013-2016.

误差估计篇7

考虑固定设计下的非参数回归模型:Yi=g (ti) +εi;i=1, 2, …n (1)

式 (1) 中A是R中的一个紧集, 固定设计点列t1, t2, …tn∈A, g (·) 是A上的有界实值未知函数, {εi;i=1, 2, …, n}为随机误差序列, 且E (εi) =0, E (ε2i) =σ2<∞, i=1, 2, …, n。不妨设A=[0, 1], 0≤t1≤t2≤…≤tn≤1。定义回归函数g (·) 的小波估计为

$\hat{g}_{n} (t) = \sum_{i = 1}^{n} Y_{i} \int_{A_{i}} E_{m} (t, s) d s (2)$

式 (2) 中Ai=[si-1, si) ;i=1, 2, …, n为区间[0, 1]上的分割且满足ti∈Ai。Em (t, s) 是由刻度函数φ (x) 产生的小波再生核:

$\begin{array}{l} E_{m} (t, s) = 2^{m} E_{0} (2^{m} t, 2^{m} s) ； \\ E_{0} (t, s) = \sum_{j \in Ζ} φ (t - j) φ (s - j) 。 \end{array}$

其中m=m (n) >0为仅依赖于n的常数。

对于模型 (1) 的回归函数的小波估计, 许多学者进行了大量的研究, 取得了丰硕的成果。如文献[1]在独立随机误差下研究了小波估计的相合性, 渐近正态性, 渐近方差;文献[2]在误差为φ混合情形下讨论了小波估计的收敛速度问题;文献[3,4]分别在α混合和ρ混合误差下研究了小波估计的相合性及收敛速度;文献[5]研究了φ混合误差下小波估计的渐进正态性;文献[6]讨论了随机误差为鞅差序列和Lq混合平稳序列时回归函数小波估计的大样本性质。

本文研究了在随机误差为 $\tilde{ρ}$ 混合时回归函数g (x) 的小波估计 (2) 的渐近正态性。

1 引理及假定

基本假定条件

(A1) 刻度函数φ (·) 是τ正则且具有紧支撑, 满足1阶Lipschitz条件, 并有

$| ϕ^{*} (ζ) - 1 | = Ο (ζ) ‚ ζ \to \infty 。$

其中ϕ*为ϕ的Fourier变换。

$(A 2) g (\cdot) \in Η^{v} ‚ v > \frac{3}{2}$ , 且g (·) 满足1阶Lipschitz条件。

$(A 3) \max_{1 \leq i \leq n} | s_{i} - s_{i - 1} | = Ο (n^{- 1})$ 。

(A4) (i) 2m=O (n1/3) , (ii) 22m/n→0。

$(A 5) \frac{2^{3 m}}{n} \to \infty$ 。

(A6) 存在正整数p:=p (n) , q:=q (n) , 使对充分大的n, 有p+q≤n, qp-1≤C<∞, 且当n→∞时

$(i) \frac{q}{p + q} 2^{m} \to 0$ ; $(i i) p \frac{2^{m}}{n} \to 0$ ;

$(i i i) k p \tilde{ρ} (q)^{1 / 2} (2^{m} / n)^{1 / 2} \to 0$ 。

假定条件 (A1) — (A4) 是讨论小波估计的一般性条件, 具体见文献[1,2,3,4,5,6]。条件 (A5) , (A6) 满足是容易验证的, 见文献[5]。

引理1 当 (A1) — (A3) 成立时, 有

(i) $\sup_{t} \int_{0}^{1}$ Em (t, s) ds<C;

(ii) $| \int_{A_{i}} E_{m} (t, s) d s | = Ο (2^{m} / n)$ ;

$\sum_{i = 1}^{n} | \int_{A_{i}} E_{m} (t, s) d s | \leq C$ ;

$\sum_{i = 1}^{n} (\int_{A_{i}} E_{m} (t, s) d s)^{2} = Ο (\frac{2^{m}}{n})$ ;

(iii) ∫ $_{0}^{1}$ Em (t, s) g (s) ds=g (t) +O (2-m) 。

(i) , (ii) 和 (iii) 的证明见文献[1,3,5,7]。

引理2[8] 当条件 (A1) — (A4) 成立时, 在模型 (1) 中若随机误差列{εi;i=1, 2…, n}为 $\tilde{ρ}$ 混合序列, 且Eεi=0;i=1, 2, …, n, 则

$E \hat{g}_{n} (t) - g (t) = Ο (2^{- m}) + Ο (n^{- 1})$ 。

引理3[8] 当条件 (A1) — (A4) 成立时, 在模型 (1) 中若随机误差列{εi;i=1, 2…, n}为 $\tilde{ρ}$ 混合平稳序列, 且 $\sum_{k = 1}^{\infty} \tilde{ρ} (k) < \infty$ , 则

$V a r (\hat{g}_{n} (t)) = Ο (\frac{2^{m}}{n})$

引理4[9] 设{Xi;i∈N}为 $\tilde{ρ}$ 混合序列, $E X_{i} = 0 ‚ E | X_{i} |^{q} < \infty ‚ q \geq 2 ‚ \tilde{ρ} (1) < 1$ , 记 $S_{n} = \sum_{i = 1}^{n} X_{i}$ , 则存在仅依赖于 $\tilde{ρ}$ (·) 和q的正常数C, 使∀n≥1有

$E | S_{n} |^{q} \leq C {\sum_{i = 1}^{n} E | X_{i} |^{q} + (\sum_{i = 1}^{n} E (X_{i})^{2})^{q / 2}}$ 。

引理5[10] 设{Xi;i∈N}为 $\tilde{ρ}$ 混合序列, p, q为两个正整数, 记

$η_{l} = \sum_{j = (l - 1) (p + q) + 1}^{(l - 1) (p + q) + p} X_{j}$ ; (1≤l≤k) 。

则有

$\begin{array}{l} | E e x p (i t \sum_{l = 1}^{k} η_{l}) - \prod_{l = 1}^{k} E \exp (i t η_{l}) | \leq \\ C | t | \tilde{ρ} (q)^{1 / 2} \sum_{l = 1}^{k} | η_{l} |_{2} 。 \end{array}$

2 主要结果及证明

定理1 当假定条件 (A1) — (A6) 成立时, 在模型 (1) 中若随机误差列{εi;i=1, 2, …, n}为 $\tilde{ρ}$ 混合同分布序列, 且 $E | ε_{i} |^{2 + δ} < \infty ‚ 0 < δ \leq 1 ‚ \tilde{ρ} (k) = Ο (n^{- θ}) ‚ θ > 1$ 则

$σ_{n}^{- 1} {\hat{g}_{n} (t) - E \hat{g}_{n} (t)} \overset{d}{\to} Ν (0, 1) ‚ t \in [0, 1] ‚ n \to \infty$ 。

证明采用文献[11]中定理2.1的证明方法, 记

$S_{n} = σ_{n}^{- 1} {\hat{g}_{n} (t) - E \hat{g}_{n} (t)} ‚$

Zni=σ $_{n}^{- 1}$ εi∫AiEm (t, s) ds, i=1, 2, …, n。

则 $S_{n} = \sum_{i = 1}^{n} Ζ_{n i}$ 。

令k=[n/ (p+q) ], 利用Bernstein大小分块原理, Sn可分解为Sn=S′n+S′n+Sn。

$\begin{array}{l} 其中 S^{'}_{n} = \sum_{m = 1}^{k} y_{n m} ‚ S^{″}_{n} = \sum_{m = 1}^{k} y^{'}_{n m} ‚ S^{‴}_{n} = y^{'}_{n k + 1} ‚ \\ y_{n m} = \sum_{i = k_{m}}^{k_{m} + p - 1} Ζ_{n i} ‚ y^{'}_{n m} = \sum_{i = l_{m}}^{l_{m} + q - 1} Ζ_{n i} ‚ y^{'}_{n k + 1} = \sum_{i = k (p + q) + 1}^{n} Ζ_{n i} ‚ \\ k_{m} = (m - 1) (p + q) + 1 ‚ \\ l_{m} = (m - 1) (p + q) + p + 1 ； m = 1, 2, \dots, k 。 \end{array}$

则在定理1的条件下, 可证得

$\begin{array}{l} S^{″}_{n} + S^{‴}_{n} \overset{Ρ}{\to} 0 (3) \\ S^{'}_{n} \overset{d}{\to} Ν (0, 1) (4) \end{array}$

成立。由式 (3) , 式 (4) 和Slutsky引理知, 定理1成立。

首先证明式 (3) 成立。由引理1 (ii) 和引理3知, $σ_{n}^{- 2} | \int_{A_{i}} E_{m} (t, s) d s | \leq C$ 。

再由引理1 (ii) 和条件 (A6) (i) :

$\begin{array}{l} E (S^{″}_{n})^{2} = E (\sum_{m = 1}^{k} \sum_{i = l_{m}}^{l_{m} + q - 1} σ_{n}^{- 1} ε_{i} \int_{A_{i}} E_{m} (t, s) d s)^{2} \leq \\ C \sum_{m = 1}^{k} \sum_{i = l_{m}}^{l_{m} + q - 1} | \int_{A_{i}} E_{m} (t, s) d s | \leq \\ C k q \frac{2^{m}}{n} \leq C \frac{n}{p + q} q \frac{2^{m}}{n} = \\ C \frac{q}{p + q} 2^{m} \to 0 (5) \end{array}$

$\begin{array}{l} E (S^{‴}_{n})^{2} = E (\sum_{i = k (p + q) + 1}^{n} σ_{n}^{- 1} ε_{i} \int_{A_{i}} E_{m} (t, s) d s)^{2} \leq \\ C \sum_{i = k (p + q) + 1}^{n} | \int_{A_{i}} E_{m} (t, s) d s | \leq C [n - k (p + q)] \frac{2^{m}}{n} \leq C (p + q) \frac{2^{m}}{n} = \\ C (1 + q p^{- 1}) p \frac{2^{m}}{n} \to 0 (6) \end{array}$

由式 (5) , 式 (6) 和Tchebychev不等式, 对∀ε>0, 有

$\begin{array}{l} Ρ (| S^{″}_{n} + S^{‴}_{n} | > 2 ε) \leq Ρ (| S^{″}_{n} | > ε) + \\ Ρ (| S^{‴}_{n} | > ε) \leq \frac{E (S^{″}_{n})^{2}}{ε^{2}} + \frac{E (S^{‴}_{n})^{2}}{ε^{2}} \to 0 。 \end{array}$

由ε的任意性知, 式 (3) 成立。

现证明式 (4) 成立。令 $s_{n}^{2} = \sum_{m = 1}^{k} V a r (y_{n m}) ‚ Γ_{n} = \sum_{1 \leq i < j \leq k} cov (y_{n i}, y_{n j})$ , 则

s $_{n}^{2}$ =E (S′n) 2-2Γn, E (S $_{n}^{2}$ ) =1,

E (S′n) 2=E[Sn- (S″n+Sn) ]2=

1+E (S″n+Sn) 2-2E[Sn (S″n+Sn) ]。

则

$| E (S^{'}_{n})^{2} - 1 | =$ E (S″n+Sn) 2-

2E[Sn (S″n+Sn) ]→0 (7)

$\begin{array}{l} | Γ_{n} | \leq \sum_{1 \leq i < j \leq k} \sum_{u = k_{i}}^{k_{i} + p - 1} \sum_{v = k_{j}}^{k_{j} + p - 1} | cov (Ζ_{n u}, Ζ_{n v}) | \leq \\ \sum_{1 \leq i < j \leq k} \sum_{u = k_{i}}^{k_{i} + p - 1} \sum_{v = k_{j}}^{k_{j} + p - 1} σ_{n}^{- 2} \end{array}$

∫AuEm (t, s) ds∫AvEm (t, s) ds $cov (ε_{u}, ε_{v}) \leq C \sum_{1 \leq i < j \leq k} \sum_{u = k_{i}}^{k_{i} + p - 1} \sum_{v = k_{j}}^{k_{j} + p - 1}$

$\begin{array}{l} | \int_{A_{u}} E_{m} (t, s) d s | \tilde{ρ} (v - u) | ε_{u} |_{2} | ε_{v} |_{2} \leq \\ C \sum_{i = 1}^{k - 1} \sum_{u = k_{i}}^{k_{i} + p - 1} | \int_{A_{u}} E_{m} (t, s) d s | \sum_{j = i + 1}^{k} \sum_{v = k_{j}}^{k_{j} + p - 1} \tilde{ρ} (v - u) \leq \end{array}$

$C \sum_{u = 1}^{n}$ ∫AuEm (t, s) ds

$\begin{array}{l} \sum_{j = q}^{\infty} \tilde{ρ} (j) \leq \\ C \sum_{j = q}^{\infty} \tilde{ρ} (j) \to 0 (q \to \infty) (8) \end{array}$

由式 (7) , 式 (8) 知

E (S′n) 2→1, s $_{n}^{2}$ →1 (9)

为了建立S′n的渐近正态性, 假设{ηnm;m=1, 2, …, k}是独立随机变量序列, 且ηnm和ynm (m=1, 2, …, k) 有相同的分布, 则有

Eηnm=0, Var (ηnm) =Var (ynm) 。

令Tnm=ηnm/sn, m=1, 2, …, k, 则

{Tnm;m=1, 2, …, k}是独立的, 且

$E Τ_{n m} = 0 ‚ \sum_{m = 1}^{k} V a r (Τ_{n m}) = 1$ 。

用ϕX (t) 表示随机变量X的特征函数, 则有

$ϕ_{\sum_{m = 1}^{k} y_{n m}} (x) - e^{- \frac{t^{2}}{2}}$ ≤

$E \exp (i t \sum_{m = 1}^{k} y_{n m}) - \prod_{m = 1}^{k} E \exp (i t y_{n m})$ +

$\prod_{m = 1}^{k} E \exp (i t y_{n m}) - e^{- \frac{t^{2}}{2}}$ ≤

$E \exp (i t \sum_{m = 1}^{k} y_{n m}) - \prod_{m = 1}^{k} E \exp (i t y_{n m})$ +

$\prod_{m = 1}^{k} E \exp (i t Τ_{n m}) - e^{- \frac{t^{2}}{2}}$ :=

I1+I2 (10)

由引理5, 引理1 (ii) 和条件 (A6) (iii) 得

$\begin{array}{l} Ι_{1} \leq C | t | \tilde{ρ} (q)^{1 / 2} \sum_{m = 1}^{k} | y_{n m} |_{2} \leq \\ C | t | \tilde{ρ} (q)^{1 / 2} \sum_{m = 1}^{k} E (\sum_{i = k_{m}}^{k_{m} + p - 1} | σ_{n}^{- 1} ε_{i} \int_{A_{i}} E_{m} (t, s) d s |^{2})^{1 / 2} \leq \\ C | t | \tilde{ρ} (q)^{1 / 2} \sum_{m = 1}^{k} \sum_{i = k_{m}}^{k_{m} + p - 1} | \int_{A_{i}} E_{m} (t, s) d s |^{1 / 2} \leq \\ C | t | \tilde{ρ} (q)^{1 / 2} k p (\frac{2^{m}}{n})^{\frac{1}{2}} \to 0 (11) \end{array}$

I2→0显然, 由此及式 (10) , 式 (11) , 可把S′n看做是独立不同分布随机变量之和。由Lyapunov中心极限定理, 要证式 (4) 成立, 只需证明存在某δ>0, 有

$\frac{1}{s_{n}^{2 + δ}} \sum_{m = 1}^{k} E | y_{n m} |^{2 + δ} \to 0 ‚ n \to \infty (12)$

由引理1, 引理3, 引理4和条件 (A6)

$\begin{array}{l} \sum_{m = 1}^{k} E | y_{n m} |^{2 + δ} \leq \\ C \sum_{m = 1}^{k} [\sum_{i = k_{m}}^{k_{m} + p - 1} E | Ζ_{n i} |^{2 + δ} + (\sum_{i = k_{m}}^{k_{m} + p - 1} E | Ζ_{n i} |^{2})^{\frac{2 + δ}{2}}] \leq \end{array}$

$\begin{array}{l} C \sum_{m = 1}^{k} {\sum_{i = k_{m}}^{k_{m} + p - 1} [σ_{n}^{- 2} | \int_{A_{i}} E_{m} (t, s) d s |^{2}]^{\frac{2 + δ}{2}} + \\ [\sum_{i = k_{m}}^{k_{m} + p - 1} | \int_{A_{i}} E_{m} (t, s) d s |]^{\frac{2 + δ}{2}}} \leq \end{array}$

$\begin{array}{l} C \sum_{m = 1}^{k} {\sum_{i = k_{m}}^{k_{m} + p - 1} | \int_{A_{i}} E_{m} (t, s) d s |^{\frac{2 + δ}{2}} + \\ [\sum_{i = k_{m}}^{k_{m} + p - 1} | \int_{A_{i}} E_{m} (t, s) d s |]^{\frac{2 + δ}{2}}} \leq \end{array}$

$\begin{array}{l} C \sum_{m = 1}^{k} {p^{\frac{δ}{2}} \sum_{i = k_{m}}^{k_{m} + p - 1} | \int_{A_{i}} E_{m} (t, s) d s |^{\frac{2 + δ}{2}} + \\ [\sum_{i = k_{m}}^{k_{m} + p - 1} | \int_{A_{i}} E_{m} (t, s) d s |]^{\frac{2 + δ}{2}}} \leq \end{array}$

$\begin{array}{l} C (p \frac{2^{m}}{n})^{\frac{δ}{2}} \sum_{i = 1}^{n} | \int_{A_{i}} E_{m} (t, s) d s | \leq \\ C (p \frac{2^{m}}{n})^{\frac{δ}{2}} \to 0 (13) \end{array}$

由式 (9) , 式 (13) 知式 (12) 成立。从而式 (4) 成立。证毕。

定理2 当假定条件 (A1) — (A6) 成立时, 且{εi}满足定理1中的所有条件, 则

$σ_{n}^{- 1} {\hat{g}_{n} (t) - g (t)} \overset{d}{\to} Ν (0, 1) ‚$

t∈[0, 1], n→∞。

$\begin{array}{l} 证明 σ_{n}^{- 1} {\hat{g}_{n} (t) - g (t)} = \\ σ_{n}^{- 1} {\hat{g}_{n} (t) - E \hat{g}_{n} (t)} + \\ σ_{n}^{- 1} {E \hat{g}_{n} (t) - g (t)} (14) \end{array}$

由定理1知 $σ_{n}^{- 1} {\hat{g}_{n} (t) - E \hat{g}_{n} (t)} \overset{d}{\to} Ν (0, 1)$ , 故只需要证明

$σ_{n}^{- 1} {E \hat{g}_{n} (t) - g (t)} \to 0 (15)$

即可。由引理2有

$E \hat{g}_{n} (t) - g (t) = Ο (2^{- m}) + Ο (n^{- 1})$ 。

再由引理3及条件 $\frac{2^{3 m}}{n} \to \infty$ 得

$\begin{array}{l} σ_{n}^{- 1} {E \hat{g}_{n} (t) - g (t)} = Ο (\sqrt{\frac{n}{2^{3 m}}}) + \\ Ο (\frac{1}{\sqrt{n 2^{m}}}) \to 0 。 \end{array}$

定理证毕。

参考文献

[1]Antoniads A, Gregoire G, Mckeague I W.Wavelet method for curve estimation.JASA, 1994;89:1340—1352

[2]薛留根.混合误差下回归函数小波估计的一致收敛速度.数学物理学报, 2002;22A (4) :528—535

[3]孙燕, 柴根象.固定设计下回归函数的小波估计.数学物理学报, 2004;24A (5) :579—606

[4]Li Yongming, Wu Lisha.Consistency of wavelet estimator of regres-sion function underρ-mixing assumptions.Journal of Shangrao Normal College, 2008;28 (6) :9—14

[5]李永明, 尹长明, 韦程东.φ混合误差下回归函数小波估计的渐进正态性.应用数学学报, 2008;31 (6) :1046—1055

[6]王二红.相依样本下非参数回归函数的小波估计.合肥:合肥工业大学, 2007

[7]Walter G G Wavelets and other orthogonal systems with applications.Florida:CRC Press, 1994

[8]潘丽静, 郭鹏江.ρ～混合误差下回归函数的小波估计.科学技术与工程, 2010;18 (10) :4363—4365

[9]杨善朝.一类随机变量部分和的矩不等式及其应用.科学通报, 1998;43 (17) :1823—1827

[10]邢国栋.ρ～混合样本下回归权函数估计的一致渐近正态性.桂林:广西师范大学, 2006

误差估计篇8

最小二乘算法是系统辨识中用得最广泛的估计方法之一。标准的递推最小二乘算法是通过极小化关于输入输出数据的一个二次准则函数,即极小化估计残差平方和而得到的算法[1]。本文基于卡尔曼滤波原理[2],从一个新的角度来推导递推最小二乘算法:定义一个参数估汁误差协方差矩阵,通过极小化该协方差矩阵而推导出递推最小二乘算法。

1最小二乘辨识算法

考虑下列单输入单输出(SISO)系统的参数估计问题,

A(z)y(t)=B(z)u(t)+v(t)(1)

(1)式中{u(t)}和{y(t)}分别是模型的输入和输出序列,{v(t)}为零均值、方差为σ2的不相关随机白噪声序列A(z)和B(z)均为单位后移算z-1的多项式[z-1y(t)=y(t-1)],且

A(z)=1+a1z-1+a2z-2+…+anaz-n,

B(z)=b1z-1+b2z-2+…+bnbz-nb。

定义参数向量θ和信息向量φ(t)分别为

θ:=[a1,a2,…,ana,b1,b2,…,bnb]T∈Rn,

φ(t):=[-y(t-1),-y(t-2),…,-y(t-na),u(t-1),u(t-2),…,u(t-nb)]T∈Rn。

那么式(1)写可成下列最小二乘辨识模型,

y(t)=φT(t)θ+v(t)(2)

本文的目标是,利用系统的输入输出数据{u(t),y(t)}或{y(t),φ(t)},基于卡尔曼滤波原理,极小化参数估计误差协方差阵,而推导辨识参数向量θ的递推最小二乘算法。

很多文献给出了计算θ估计白 $\hat{θ}$ (t)的递推最小二乘算法[1,3]:

$\hat{θ} (t) = \hat{θ} (t - 1) + L (t) [y (t) - φ^{Τ} (t) \hat{θ} (t - 1)]$ (3)

L(t)=P(t-1)φ(t)[1+φT(t)P(t-1)φ(t)]-1(4)

P(t)=[I-L(t)φT(t)]P(t-1)(5)

其中 $\hat{θ}$ (t)为t时刻的估计,L(t)称为增益向量,P(t)称为协方差矩阵。

2算法推导

2.1单输入单输出系统的递推最小二乘算法

仍考虑上述单输入单输出系统式(1)或式(2)。参照式(3),设t寸刻的参数估计 $\hat{θ}$ (t)足由前一时刻参数估计 $\hat{θ} (t - 1)$ 加上一个补偿项得到,且具有下列形,

$\hat{θ} (t) = \hat{θ} (t - 1) + L (t) [y (t) - φ^{Τ} (t) \hat{θ} (t - 1)]$ (6)

那么问题就是如何找到增益向量L(t)∈Rn。方法如下,定义参数估计误差为

$\hat{θ}$ (t): $= \hat{θ} (t) - θ$ (7)

将式(6)代入式(7),并利用式(2)展开可得

$\hat{θ} (t) = [Ι - L (t) φ^{Τ} (t)] \hat{θ} (t - 1) + L (t) v (t)$ (8)

(8)式中I是一个适当维数的单位阵.目标是确定一个最优增益向量L(t)使参数估计误差 $\hat{θ}$ (t)最小。由于v(t)是白噪声,假设L(t)和φ(t)与v(t)独立,对式(8)两边取期望得 $E [\hat{θ} (t)] = 0 (E$ 为期望算子)。定义参数估计误差协方差阵.

P(t): $= E [\hat{θ} (t) \hat{θ}^{Τ} (t)]$ (9)

将式(8)代入式(9)可得

P(t)=[I-L(t)φT(t)]P(t-1)×

[I-L(t)φT(t)]T+L(t)σ2LT(t)(10)

因为P(t)为非负定矩阵,将P(t)配成下列形式:

把P(t)看作待定增益向量L(t)的函数。通过极小化估计误差协方差矩P(t),可求得最优增益向量L(t).上式中矩阵P(t)包含3项,第一、二项与L(t)无关,第二项中,因为P(t-1)是非负定矩阵,所以σ2+φT(t)P(t-1)φ(t)≥0。如果选择增益L(t)使得第三项为零,即取

L(t)=P(t-1)φ(t)[σ2+φT(t)P(t-1)φ(t)]-1(12)

那么协方差阵P(t)最小,也就是参数估计误差最小。

$E [∥ \hat{θ} (t) ∥^{2}] = t r [Ρ (t)] = \min$ 。

将式(12)代入式(11)得到

P(t)=P(t-1)-P(t-1)φ(t)×

[σ2+φT(t)P(t-1)φ(t)]-1φT(t)P(t-1)=

[I-L(t)φT(t)]P(t-1)(13)

式(6),式(12)和式(13)构成了基于卡尔曼原理的单输入单输出系统的最小方差递推最小二乘算法它与标准递推最小二乘算法式(3)-式(5)的差别在于式(12)-式(13)中取σ2=1.

下面讨论一种特殊情况。对于确定性系统,噪声方差σ2=0,则算法式(6),式(12)-式(13)退化为

$\hat{θ} (t) = \hat{θ} (t - 1) + L (t) [y (t) - φ^{Τ} (t) \hat{θ} (t - 1)]$ ,(14)

L(t)=P(t-1)φ(t)[φT(t)P(t-1)φ(t)]-1(15)

P(t)=[I-L(t)φT(t)]P(t-1)(16)

为了防止式(15)中分母为零,可加上一个小常数ε(如取ε=10-6),则算法式(14)-式(16)可修改为

$\hat{θ} (t) = \hat{θ} (t - 1) + L (t) [y (t) - φ^{Τ} (t) \hat{θ} (t - 1)]$ (17)

$L (t) = \frac{Ρ (t - 1) φ (t)}{ε + φ^{Τ} (t) Ρ (t - 1) φ (t)}, 0 ＜ ε \leq 1$ (18)

P(t)=[I-L(t)φT(t)]P(t-1)(19)

在后面的例子中,将研究ε取不同值时,参数估计精度的变化。

这是文献[2]状态估计算法推导过程,推导递推最小二乘参数估计的一种新的简便方法。从上述推导看:这种参数状态估计误差具有最小方差性质。也可理论上证明这个递推最小二来算法的参数估计误差收敛于零。

2.2多输入多输出系统的递推最小二乘算法

上述方法也叮推广到多输入多输出(MIMO)系统.考虑下列多输入多输出系统,

A(z)y(t)=B(z)u(t)+v(t)(20)

A(z):=I+A1z-1+A2z-2+…+Anaz-na,

B(z):=B1z-1+B2z-2+…+Bnbz-nb,

u(t)∈Rr为系统输入向量,y(t)∈Rm为系统输出向量,u(t)∈Rm为零均值、协方差阵为Rv=E[v(t)vT(t)]∈Rm×m的随机噪声向量。

定义参数矩阵ϑ和信号向量φ(t)如下

ϑT:=[A1,A2,…,Ana,B1,…,Bnb]∈Rm×m,φ(t):=[-yT(t-1),…,-yT(t-ma),tT(t-1),…,uT(t-nb)]T∈Rn。

则式(20)可等价写为

y(t)=ϑTϕ(t)+v(t)(21)

参照式(6),由上式可知ϑ的递推估计应取为下列形式:

$\hat{ϑ} (t) = \hat{ϑ} (t - 1) + L (t) [y^{Τ} (t) - ϕ^{Τ} (t) \hat{ϑ} (t - 1)]$ (22)

则参数估计误差可写为

$\begin{array}{l} \hat{ϑ} (t) ∶ = \hat{ϑ} (t) - ϑ = [Ι - L (t) ϕ^{Τ} (t)] \times \\ \hat{ϑ} (t - 1) + L (t) v^{Τ} (t) . \end{array}$

定义参数估计误差协方差短阵P(t): $= E [\hat{ϑ} (t) \times \hat{ϑ}^{Τ} (t)]$ 。令ϵ:=tr[Rv]=E[‖v(t)‖2].采用与上节相同的方法,可以得到多输入多输出系统的最小方差递推最小二乘算法:

$\hat{ϑ} (t) = \hat{ϑ} (t - 1) + L (t) [y^{Τ} (t) - ϕ^{Τ} (t) \hat{ϑ} (t - 1)]$ (23)

$L (t) = \frac{Ρ (t - 1) ϕ (t)}{ϵ + ϕ^{Τ} (t) Ρ (t - 1) ϕ (t)}$ (24)

P(t)=[I-L(t)ϕT(t)]P(t-1)](25)

如果式(24)中取ϵ=1,就得到多输入多输出系统的标准递推最小二乘算法.

3仿真例子

例1 考虑下列随机系统,

A(z)y(t)=B(z)u(t)+v(t),

A(z)=1+a1z-1+a2z-2=1-1.60z-1+0.80z-2,

B(z)=b1z-1+b2z-2=0.412z-1+0.309z-2,

θ=[a1,a2,b1,b2]T=[-1.60,0.80,0.412,0.309]T。

仿真时,输入u(t)采用零均值单位方差不相关可测随机序列,v(t)采用零均值方差为σ2的白噪声序列.当噪声方差为σ2=0.102和σ2=1.002时,对应的噪信比分别为δns=14.26%和δns=142.63%.应用提出的算法式(6),式(12)-式(13)估计这个系统参数,不同噪声方差下的参数估计及其误差如表1所示,参数估计误差δ: $= ∥ \hat{θ} (t) - θ ∥ / ∥ θ ∥$ 随t变化曲线如图1所示.