第2讲：捷联惯导系统(1-2)时间序列建模

3．随机漂移误差的时间序列分析建模在惯性导航系统中，为了减小陀螺随机误差对系统精度的影响，有效可行的办法是采用滤波技术对随机误差进行实时补偿。实时补偿的前提是已知随机误差的模型，为此，需要事先对陀螺的随机噪声进行必要的数学处理，建立适合于在线补偿的数学模型。陀螺的随机噪声一般是有色噪声，即非平稳的随机过程，处理这类随机过程较成熟的建模方法是时间序列分析法。该方法是针对一组离散随机数据序列，进行时域和频域内的统计特性分析，求出实际物理系统的统计特性，并将随机数据浓缩成一个简单的随机差分模型。时间序列分析建模流程图3.1序列的检验随机数据处理分析的结果是否正确，取决于数据的一些基本特性。其主要的三个基本特性是：平稳性；周期性；正态性。对这三个基本特性进行检验，是分析、建立光学陀螺随机误差模型的重要前提。1．平稳性检验主要目的是检验陀螺随机误差时间序列是否具有不随时间原点的推移而变化的统计特性。如果陀螺随机误差时间序列是平稳的，再加上假设为各态历经的，则对陀螺随机误差的研究，就可以用单个样本记录的时间序列来代替总体平均。这就给数据处理带来了极大的方便。如果不平稳，则需要对数据进行平稳化处理。造成随机过程不平稳的原因，是随机过程中包含有随时间缓慢变化的趋势项。检验这种非平稳趋势项的一种很有效的方法是逆序法。（1）逆序法对于测试数据记录，将其分成nyyy,,,21LM段，然后求各段的均值（或方差值），得到一个大致不相关的均值（或方差值）序列Mxxx,,,21L对于下标为的，每当出现iixijxx)1,,2,1,(−=MiijL时就定义为的一个逆序，与相应的逆序的个数称为的逆序数。序列ixixiAix121,,,−MxxxL的逆序总数定义为∑−==11MiiAA以随机整数序列出现的的均值与方差分别为A[][]()41211111−===∑∑−=−=MMiAEAEMiMii[]()725322−+=MMMAVar统计量[][]AVarAEAu⎟⎠⎞⎜⎝⎛−+=21渐近服从正态分布)1,0(N。平稳性的判断：如果值处在之内，则可接受序列无趋势的假设（在0.05显著水平上）；否则拒绝该假设。u2±A很大时，表明序列均值（或方差）有上升趋势，很小时，表明序列均值（或方差）有下降趋势。A（2）平稳化处理对随机序列进行平稳性检验，如发现为非平稳的随机序列，则需要对其进行平稳化处理，即提取趋势项。趋势项的含义是数据的周期大于样本长度的项。常规处理方法是用线性项去拟合趋势项。但实际上，趋势项并不一定是标准的直线，而有可能是一个随均值变化的曲线。因此，一般将趋势项表示为时间t的多项式()nntatataatf++++=.....2210然后用线性回归的方法，拟合各项系数。考虑非平稳随机序列中可能存在多种长期变化趋势，也可用如下函数来近似拟合()tettttttttftln1092182172615432210αααααααααα++++++++=−−−−经过趋势项提取处理后得到的残差序列，仍有可能是不平稳的时间序列（尤其是对于光学陀螺）。这主要是由于随机游走造成的，必须用差分的方法来处理。对于含有随机游走成分的非平稳时间序列，只需经过一阶差分，即可化为平稳时间序列。在对非平稳时间序列进行差分处理时，必须注意差分处理会使系统的频率特性发生畸变。因为微分环节将使系统的低频部分发生严重的衰减。因此，在对陀螺随机误差数据序列进行平稳化处理时应尽量少用差分处理，差分次数不应超过两次。2．周期性检验周期性检验就是用来识别陀螺随机误差数据中是否含有随机量以外的周期性分量，周期性检验的方法是直接考察从陀螺随机误差数据得到的自相关函数或功率谱密度的图形。（1）相关函数识别随机序列中如果存在周期项，会在自相关函数图形中反映出来。随机量和周期量的自相关函数图形有明显的不同：a.随机量的图形在时间间隔τ增大时，总是一条衰减的曲线；b.周期量的图形不管时间间隔τ怎样增大，总是一条不衰减的振荡曲线；c.如果随机序列中同时含有随机量和周期量，就会出现自相关函数曲线在一定的时间间隔τ内呈衰减趋势，然后便维持不衰减的振荡。（2）功率谱密度识别当数据中同时含有随机量和周期量时，就会出现功率谱密度曲线中有明显的尖峰。但是如果随机序列中周期性不是太强，则其尖峰可能不很明显。（3）周期图分析法前述两种方法只能对数据中是否含有周期项进行判断。为了寻找隐含周期，识别并提取随机序列中的周期函数项，通常采用周期图分析法。a.周期图的定义假定平稳的时间序列Nnxxxx,,,,21LL是由各种不同频率的正弦波、余弦波叠加而成的函数，对其进行傅立叶展开，有NtNjtbNjtaaxkjjjt,,2,12sin2cos10L=⎟⎠⎞⎜⎝⎛++=∑=ππ式中⎪⎩⎪⎨⎧−=为奇数时当为偶数时当NNNNk212傅立叶系数为⎟⎠⎞⎜⎝⎛⎥⎦⎤⎢⎣⎡−=====∑∑∑===21,,2,12sin22cos211110NjNjtxNbNjtxNaxxNaNttjNttjNttLππ定义()⎥⎦⎤⎢⎣⎡−=+=21,,2,12)(22NjbaNfIjjjL为时间序列的周期图。b.隐含周期的提取对于随机序列Nnxxxx,,,,21LL，考虑含有隐含周期项的模型为)()sincos()(10ttbtaatxliiiiiηωω+++=∑=其中，),,2,1(,,,0libaaiiiL=ω和都是有待识别、估计的参数。l计算周期图)(jfI),,2,1(KjL=。下面利用Fisher统计检验方法从周期jNTj=中识别出时间序列的真正可能周期。（1）第1个周期的识别与提取取{})()(,),(),(max1121jKjfIfIfIfII==L在)(tη为正态白噪声的假定下，统计量∑==KjjjfIIg1)1(1服从Fisher分布{}[]∑=−++−−=rjKjKjgjCggp01111)1(1)1(其中，是使成立的昀大正整数。r0)1(11+−gr给定显著水平α：若{}α≥1ggp，可以断言时间序列无隐含周期函数项；若{}α1ggp，可以接受11jNT=为时间序列的第一个周期，相应的傅立叶系数给出了隐含周期模型中的系数的估计值。11jjba、11ba、（2）第k个周期的识别与提取取kjI为周期图)(,),(),(21KfIfIfIL中的第()1kk个昀大值∑==KjjjkfIfIgk1)()(Fisher检验中的统计量定义为∑−=−=′111kjjkkggg对给定的显著水平α：若{}α′kggp成立，则接受kkjNT=为时间序列的一个周期。显著水平α不宜取得过大，一般取01.0≤α，否则可能选入较多的伪周期。3．正态性检验如果随机误差时间序列符合平稳性检验，又符合正态性检验，则该随机过程必定是满足各态遍历性的平稳随机过程。在实际工作中，昀常用的正态性检验的方法是拟合优度检验法，实现思想如下。2χ(1)求观测频数对于实验数据Nxxx,,,21L，将它的一个区间范围bxa分成k个相等的子区间，整个的范围看成由x2+k个区间构成。令kabc−=，jcadj+=，(kj,,1,0L=)以iN表示出现在第区间的数据点个数，称为观测频数，则iiN的计算如下1、若，则使axi0N加1；2、若，则使bxi≥1+kN加1；3、若上述两个条件均不满足，则bxai≤，计算1+⎥⎦⎤⎢⎣⎡−=caxji，求得j后令jN加1（运算符[]•表示取实数的整数部分）。利用110,,,+kNNNL所画出的矩形阶梯状图形就称为概率直方图。(2)正态性检验——拟合优度检验2χ正态分布的概率密度函数为()22221)(σµσπ−−=xexp式中，µ、σ分别为总体的均值与方差。概率分布函数是概率密度函数的积分()()⎟⎠⎞⎜⎝⎛−≡=∫∞−−−σµφσπσµxdteXxPxt22221其中，φ称为概率积分。随机变量处于α和β之间的概率为()⎟⎠⎞⎜⎝⎛−−⎟⎠⎞⎜⎝⎛−=≤σµαφσµβφβαxP对于序列Nxxx,,,21L，如果数据分布是正态的，则可计算处应落入第j个子区间的数据个数（称为期望频数）为⎭⎬⎫⎩⎨⎧⎟⎠⎞⎜⎝⎛−−=⎭⎬⎫⎩⎨⎧⎟⎟⎠⎞⎜⎜⎝⎛−−⎟⎟⎠⎞⎜⎜⎝⎛−=⎟⎠⎞⎜⎝⎛−=+−σµφσµφσµφσµφˆˆ1ˆˆˆˆˆˆ110bNFddNFaNFkjjjLLLLLL均值与方差加上标“∧”表示样本均值与样本方差。根据观测频数jN构造统计量()∑+=−=1022kjjjjFFNχ假设这个统计量渐近服从自由度为1−k的分布，将其与理论的分布进行比较。2χ2χ设给定的显著水平为α：若，则接受数据2,12αχχ−≤kNxxx,,,21L为正态分布的序列；否则，认为该数据不是正态分布的序列。3.2平稳序列的建模经过平稳化处理以及周期项分析与提取，得到了平稳的随机漂移序列，对于平稳随机序列可用线性的AR、MA或ARMA模型进行拟合。即一个平稳的随机过程可看成是不相关的白噪声通过上述三种中的任意一种模型而产生的过程。ARMA(,)是AR模型和MA模型的混合模型，它的一般结构式为pqqtqtttptpttt−−−−−−+⋅⋅⋅+++=+⋅⋅⋅+++θθθφφφ22112211式中，pφφφ,......,,21为自相关系数，qθθθ,......,21为滑动平均系数。ARMA(,)描述的是这样一种平稳随机过程：该随机过程在时刻的观测值与时刻之前的观测值存在相关性，而且时刻的观测值与时刻之前的个白噪声也存在相关性。pqttxtptttxxx−−−,......,21ttxtqqttt−−,,,1L（1）确定模型结构对于平稳的随机序列，可利用其自相关函数和偏自相关函数的特点来初步识别模型的形式。对于平稳的随机序列Nηηη,,,21L，自相关函数的计算公式如下：Nmmkrrkk==;,,2,1,0ˆˆˆ0Lρ其中，∑−=+−−=kNtkttkNr1))((1ˆηηηη。偏自相关函数的递推计算公式：),,2,1(ˆˆ1ˆˆˆˆˆˆˆˆˆˆ11111,1)1(,1,1,1111kjkjkjjkjkjjkkkkjkkkkkjjkL=−−=−==∑∑==+−+++−−+++ϕρϕρρϕϕϕϕϕρϕ三种线性平稳模型的相关特性)(pAR)(qMA),(qpARMA自相关函数kρ拖尾截尾)(qk拖尾偏相关函数mmϕ截尾)(pm拖尾拖尾（2）模型参数估计时间序列线性模型的参数估计，是指在辩识得到模型类型和阶数的基础上，求出模型中自相关系数及滑动平均系数的值。通常采用的方法有昀小二乘法和极大似然法。（3）模型适用性检验经前面的分析与处理，可知陀螺随机漂移的拟合模型为ttttxfyηˆˆˆˆ++=其中，为趋势项，为周期项，tfˆtxˆtηˆ为自回归模型。这里加上上标‘∧’以区别于真实值。模型检验的主要依据是陀螺随机漂移拟合残差taˆtttyyaˆˆ−=taˆ愈接近白噪声说明模型的拟合精度愈高。检验拟合残差是否为白噪声序列也就是判断的独立性假设成立与否。taˆtaˆ白噪声检验方法：对于白噪声序列，当1≥k时，样本自相关函数的分布渐近于正态，即()NNr1,0~)(ˆρ所以，当N较大时{})(ˆ,),2(ˆ),1(ˆkNNNρρρL这k个量可近似为相互独立的满足正态分布()1,0N的随机变量，因而它们的平方和符合分布。构造统计量2χ∑==krrNQ12)(ˆρ检验残差序列是否为白噪声的问题，就转化为检验统计量是否是自由度为taˆQk的分布的问题。2χ3.3建模实例实测的某型陀螺上电启动时所记录的一段输出数据如下图。01002003004005000.02080.0210.02120.02140.02160.02180.0220.02220.0224测试数据测试数据曲线（1）数据的统计特性对测试数据进行统计分析，其均值与方差分别为：均值=0.0220，方差=4.1781e-008昀大的相关延迟取为5/N，则相关函数与偏相关函数计算结果如下：相关函数：1.0000，0.9696，0.9410，0.9142，0.88940.866