非参数估计法（精选九篇）
非参数估计法 篇1
资本资产定价模型是金融学的基石, 同时也是学术界研究最多, 争论最多的理论。在金融资产定价模型中, 很多都是预测资产收益模型, 如:资本套利模型、基于消费的均衡模型。但是, 没有一个模型能够像Sharpe-Lintner的条件CAPM模型一样受学术界的青睐。CAPM模型是建立在市场组合均值—方差有效的假定基础之上, 并且在这一假设下认为单个风险资产的收益与市场资产组合的风险收益是成比例, 其中β为市场有价证券的系数, 用来衡量市场有价证券收益对市场风险变动的敏感程度。这个简单的CAPM模型就是众所周知的无条件或者是静态CAPM模型, 在这个模型里, 单个有价证券和市场资产组合的关系是不随时间变化的, 也既是β不随时间和市场波动而变化。
在过去的几十年里, 学者们对CAPM模型进行了大量的实证检验, 静态CAPM模型的许多异像被发现。Fmam—French (1992) 提出静态CAPM不支持实证研究的观点, 就像重磅炸弹一样在理论界和实业界引起震动, 很多人对CAPM模型的信心开始动摇, 甚至有人认为CAPM已经死亡。但是, 仍然有很多学者是支持CAPM, 他们为此进行着不懈的努力, 有部分学者将注意力放在了β稳定性方面, Levy建议分市场研究β, Fabozzi和Francis分别对牛市和熊市的β稳定性作了检验。他们发现资产定价模型中的单个市场指数是不受牛市和熊市影响的。
另一方面, Keim和Stambaaugh, Breen, Glosten和Jagannathan认为在CAPM框架中β不是静态的, 而是时变的。Chen, Ferson和Harvey也提出了β是随商业周期而变化的。在Jagannathan和Wang的 (1996) 论文中拓展了条件CAPM模型, 在该条件CAPM模型中有价证券的β是由投资者在t时刻可利用的信息集而决定的, 并且随着经济情况的波动而变化。
条件CAPM的发展激发了学者们又把焦点放在了对条件模型的形成和检测方面。尽管条件CAPM能够对静态CAPM的异像提出一定的解决方法, 但其本身也产生了一些新的问题, 其中一个问题就是对变动因素的选择以及β与各个变动因素之间究竟是什么样的关系缺少理论的支持。最初, 有些学者以β与变动因素之间是线性的函数关系来进行实证检验。然而, 这种检验方法的结果有时会得到比静态CAPM模型更糟糕的结果。Ghysels认为条件CAPM定价错误的原因就在于人们认为β与动态风险之间的函数关系像静态CAPM模型一样是线性的函数关系导致的。
为了解决条件CAPM在实证中的问题, 很多学者把眼光放在了无参数估计技术方面, 采用非参技术可以避免采用β和变动因素原有的特定假设函数形式, 从而提高检验的准确度。王振宇在他的文章中提出了一种新的灵活的非参数检验方法, 该方法建立的基础是对隐含于条件线性因子定价模型中的随机折现因子的非参数限制。在检验中该方法脱离了对条件β, 风险升水和随机折现因子原有的函数形式。本文正是利用王振宇提出的该非参数检验方法利用中国沪市A股数据对条件CAPM模型进行实证检验, 验证中国股市是否存在公司规模和账面市值比效应, 条件CAPM模型在中国股市是否成立。
二、检验方法的理论基础
条件资本定价模型形如:
Ri, t表示均衡状态下证券i在t时刻的收益率变量, RM, t表示市场组合证券在t时刻的收益率变量, Rf为无风险收益率, It-1表是t-1时刻所有与风险资产价格相关的信息集。条件资本资产定价模型是将静态的资本资产定价模型中的风险资产收益、市场组合收益率变量增加条件限制, 假设他们的变化受前期信息集的影响, 在这种定义下β系数也就不再是固定的, 而是随前期信息或其他变量信息的变动而变动。这样, 模型对预期收益的解释程度便会随之加强。
如前所述, 王振宇的非参数检验方法是依赖于对隐含于条件资本资产定价模型中的随机折现因子框架的限制之上, 随机折现因子框架非常通用的框架。其对任何现代资本资产定价模型都成立的基底方程为:
其中Et表示条件收益, mt+1表示随机折现因子, Ri, t+1表示资产i的收益。
方程 (2) 也等价于下式:
n表示资产的个数, ri, t+1表示资产i的超额收益。
对于方程 (1)
因此有
为了实证目的, 令xt为条件变量集, 且, ,
这里超额收益、条件变量假定为严格静态的。
在 (5) 式假设下, 条件资本资产定价模型的条件定价误差为:
这里mt+1=1-b (xt) rp, t+1, 与 (3) Et (mt+1ri, t+1) =0, i=1, …N表示意思相同。本文采用与Wang相同的Nadaraya-Watson核估计方法来估计非参数的随机折现因子。
其中K (·) 为核函数, h为窗宽。
则Nadaraya-Watson核回归函数为:
摘要：资本资产定价模型是金融学的基石, 国内现有的对条件CAPM的实证研究大部分采用的是参数估计方法, 在做实证研究时。本文着重介绍非参数估计理论, 以解决条件CAPM在实证中的问题, 从而提高检验的准确度。
关键词：静态CAPM,条件CAPM,随机折现因子,核函数
参考文献
[1]Bure Kayahan, Thanasis Stengos, Testing the capital asset pricing model with Local Maximum Likelihood methods[J].Seience Direct46 (2007) 138～150
[2]Turan G.Bali, The intertemporal relation between expected returns and risk[J].Journal of Financial Economics87 (2008) 101–131
[3]曹培慎:金融资产定价的随机折现因子方法[J].统计与决策, 2007.07 (理论版)
非参数估计法 篇2
α混合下非参数回归函数改良分割估计的强相合性及收敛速度
In this paper, we study the strong consistency and convergence rate of modified partitioning estimate of nonparametric regression function under the sample {(Xi, Yi),i ≥ 1} that is α sequence taking values in Rd × R1 with identical distribution.
作 者：姚梅 杜雪樵 YAO Mei DU Xue Qiao 作者单位：Department of Mathematics, Hefei University of Technology, Anhui 230009, China 刊 名：数学研究与评论 ISTIC PKU英文刊名：JOURNAL OF MATHEMATICAL RESEARCH AND EXPOSITION 年，卷(期)： 28(3) 分类号：O212.7 关键词：nonparametric regression function modified partitioning estimate strong consistency convergence rate α-mixing基于相位法的雷达脉内参数估计 篇3
关键词：瞬时测频；多项式拟合法；脉内参数估计
雷达信号的脉内特征是雷达信号相位特征的重要体现，是电子侦察中对雷达信号分选和识别的重要参数。因此，要可靠的分选和识别雷达信号就必须对雷达信号进行脉内特征分析。
在各种脉内分析方法中，分析信号的瞬时频率变化规律是进行脉内调制分析的有效手段。由于信号的频率是信号相位的变化率，所以对信号相位的估计自然可以用来对信号瞬时频率进行估计，从而进行相关脉内参数的估计，因此该方法的前提是能够得到准确的瞬时频率。
文中以线性调频信号以及相位编码信号为例，从侦察接收机的信号相位入手，利用瞬时相位法实现对瞬时频率的快速、准确的推算，然后利用多项式拟合法[1]来估计信号的相位参数。在Matlab中进行了仿真验证，同时给出了在某电子侦察系统中工程实现的过程。
4 结束语
本文讨论了利用信号相位建模进行瞬时频率估计算法，将频率的估计转化为求解信号模型参数的问题。在满足白噪声的情况下，由于在参数的估计过程中使用了最小二乘方法，因此能够有效的消除噪声的干扰，获得对信号瞬时频率的准确估计。基于瞬时频率进行脉内信号的调制特征的识别，还能够实现对常规脉冲、线性调频信号、非线性调频信号以及相位编码信号等调制类型信号的识别，提取相应信号的调制参数。该方法在一定信噪比条件下有较高的正确识别率，对脉宽变化不敏感，更加符合雷达侦察数字接收机上高速实现，并且在某雷达侦察数字接收机系统中已经取得成功验证，但是该类算法的问题在于，必须精确的建立相位模型，否则估计性能会受到很大的影响。
参考文献
[1] 林琳. 多项式相位信号的参数估计.实验科学与技术. 2006年6月第三期.
[2] 张春杰，郜丽鹏，司锡才. 瞬时相位法线性调频信号瞬时频率提取技术研究，弹箭与制导学报 TN911 2006.
[3] 齐国清，贾欣乐. 插值FFT估计正弦信号频率的精度分析. 电子学报 Vol.32 No.4 Apr. 2004.
[4] 胡来招. 数字瞬时测频—相位推算法测频[ J] . 电子对抗，2005.
非参数估计法 篇4
一、GARCH模型的表现形式
对数收益率序列, 我们假定其均值方程是一个ARMA模型, 设αt=rt-μt是均值修正的对数收益率, 若:
其中是{εt}一个独立同分布的随机变量序列, 均值为0, 方差为1, 参数 , (这里i>m, αi=0对j>sβj=0, ) , 则称αt服从GARCH (m, s) 模型。对αi+βj的限制条件保证的无条件方差是有限的, 同时它的条件方差σ2i是随时间变化的。如前面一样, 通常假设是标准正态分布或标准化的学生-t分布。若s=0, 则上述方程就变为一个纯ARCH (m) 模型。
二、非参数GRACH (p, q) 模型思想及估计算法
设{rt}为一平稳时间序列, Ft-1=σ (rt-1, rt-2, ......) 表示由直到时刻的过去信息产生的σ-域, rt的非参数GARCH (p, q) 模型可描述为:
其中{εt:t∈Z}为独立同分布的随机误差序列, 均值为零方差为, 并且存在有限四阶矩。εi与{rs:s≤t}相互独立, 并且的严格正值函数:σ2t=Var[rt|Ft-1]就是的条件方差, σt就是通常所说的波动率。
三、非参数GARCH模型的具体应用步骤
将3改写成具有可加性误差形式的模型:
其中为一个鞅差序列, 其均值E[Vt]=vt|Ft-1]=0;时s
(1) 对一般的标准参数GARCH (p, q) 模型运用极大似然法计算波动率的第一次估计 , 令m=1。
(2) 以 为权数, 用 对 和 做作p+q个变量的非参数加权回归, 获得f的估计f。
(3) 计算 , 边界值σ21, m用σ2, m-1代替。
(4) 增加m, 如果, 回到第 (2) 步, 其中M为事先设定的最大迭代次数。
(5) 考虑最后K步波动率的平均 , 然后作εt2对rt-1, rt-2, ...rt-p和ε2t-1, ......ε2t-q的非参数回归, 获得函数的最终估计以及条件以方差的最终估计形式。
四、非参数GARCH模型和GARCH模型对比
非参数GARCH模型具有广泛的应用范围, 同时, 非参数GARCH模型放宽了参数GARCH模型的很多限制, 例如模型的形式和误差分布, 为波动性问题的研究提供了非常有力的工具。使用非参数GARCH模型来拟合数据可以减少模型由于设定误差带来的偏差, 在较大的范围内拟合模型, 寻找更适合的线性形式。GARCH模型簇已经成为度量金融市场波动性的强有力工具。许多专家和学者做了许多关于各种GARCH模型的理论及应用、波动率预测应用。可以说在研究股市波动方面, 非参数GARCH模型已被证明是非常成功的。在其他方面, 非参数GARCH模型的优势还有待发掘。尤其是我国基金市场, 其波动性和股市的波动性有非常紧密的联系, 具有许多相似的性质, 因此我们利用非参数GARCH模型研究我国基金等投资市场的波动性具有很好前景。
参考文献
[1]Engle, R.F..Autoregressive conditional heteroscedasticity with estimates of the variance of United Kingdom inflation[J].Econometrica, 1982, 50, 987-1008
[2]Bollerslev, T..Generalized autoregression conditional heteroskedasticity[J].Journal of Econometrics, 1986, 31, 307-327
[3]Xiangdong Long Aman Ullah.Nonparametric and Semiparametric Multivariate GARCH Model[J]
[4]鲁万波:基于非参数GARCH模型的中国股市波动性预测[J].数理统计与管理, 2006, 4 (25) , 455-461
[5]许冰, 任军峰:基于非参数GARCH一种波动率估计方法[J].统计与决策, 2006, 10, 6-10
非参数估计法 篇5
1972年, Z.Jelinski和P.Moranda 提出的JM模型[4]是软件可靠性模型的代表, 它是面向时间的模型, 属于出错计数模型, 既可以估计软件的可靠性, 又可以给出程序中残留错误的数目, 但是它的假设排错是完全的, 不能适应由于排错不完全引起的新错的情况, 所以, 这并不符合实际情况。1986年出现多故障模型, 它允许在同一时刻, 查出、改正或引入的错误个数可以多于一个[5]。但这些假设与处理还没有完全说明实际情况, 在改错过程当中, 引入错误个数与剩余错误个数并不一定存在正比关系, 改正不同的错误会对软件产生不同的影响。所以, 有必要针对故障率与剩余错误个数成其他关系的情形进行研究。
对故障率与剩余错误个数成指数关系的情形进行了研究, 在排错不完全的假设下对软件的初始错误个数和指数系数给出了一个合理的估计方法。
1 模型的建立以及参数估计
1.1 模型的假设
1) 软件中初始错误个数为一个未知但固定的常数, 用X表示;
2) 将错误引入数看成负的消除数, 设软件错误消除数服从参数为θ的指数分布;
3) 设μ表示第i次改正错误个数的数学期望, σ表示改正错误个数的方差。
其中:
μ>0表明测试过程是顺利进行的, 软件总错误个数在减少, 软件可靠性随之提高。
μ=0表明测试工作没有收效, 软件总错误个数基本不变, 软件可靠性基本不变。
μ<0表明测试工作是失败的, 软件总错误个数在增加, 软件可靠性随之降低。
σ值的大小表示测试工作过程中改正错误数波动的大小[7]。
μ, σ的值可以由测试记录数据估计得出。设ai是每次改正的消除数, 则
undefined
4) 错误一旦查出, 立即改正错误。分三种情况考虑, 这三种情况发生的概率相同, 服从离散型均匀分布。
情况一:发现错误就能顺利改正错误, 于是每次改正之后, X就要减去μ;
情况二:发现错误改正错误的同时又产生一个错误, 于是每次改正之后, X就要加上μ;
情况三:发现错误改正错误的同时又产生多个错误, 于是每次改正之后, X就要加上 (μ+a) 。
5) 在任何时候的故障率都与软件中的剩余错误个数成指数关系, 指数系数用θ表示。
1.2 模型的建立
根据假设, 在第一个错误被排除之后, 故障率会发生变化。
情况一:故障率由undefined变为undefined;
情况二:故障率由undefined变为undefined;
情况三:故障率由undefined变为undefined。
设t1, t2, …, tn表示相继出现的并被改正的错误之间的时间区间样本。根据假设, 在每两个错误之间只有唯一的故障率, 关于ti的密度即为:
1.3 参数的估计
模型中含有两个未知参数X和θ, 所以, 对模型的参数估计就是X和θ的值。
由 (1) 式可得到对数似然函数为:
undefined
令undefined, 对 (2) 式中的X和θ求偏导, 并令结果为0, 得
undefined
undefined (4)
方程 (3) 式中都不含θ, 而由测试收集的数据, 可计算出undefined与undefined的值。将它们代入 (3) 式中, 则可以解出X的估计值undefined:
undefined
再将X代入 (4) 式中, 可以解出θ的估计值undefined:
undefined
2 模型评价
模型 (下称优化模型) 在以下几个方面有所改进:
1) 软件所改正的错误数与软件错误数的改正的次数两个概念完全分开, 并且软件所改正的错误数大于软件错误数的改正的次数。
2) 在任何情况下, 改正错误时既可能减少错误, 也可能增加错误。
以纵轴表示改错过程中软件剩余的错误个数, 横轴表示时间, 某次测试的优化模型剩余的错误个数与时间关系, 如图1所示。
3) 从图1可以看出, 优化模型的剩余错误量都是在减少的, 但改正错误并不一定会使错误总数减少, 反而可能出现错误数增多的现象, 这与实际情况是比较相吻合的。
4) 当时, 它就是原来的JM模型。
在研究中, 仅考虑了三种情况为离散型均匀分布, 即等可能分布, 在以后的研究中, 可以进一步研究不是等可能分布的情况。
摘要：20世纪70年代提出的JM模型是最早代表软件可靠性的数学模型, 其假设排错是完全的, 这与实际情况不相符合。因此, 有必要对JM模型假设进行改进, 从而提高软件可靠性。在排错是不完全的假设下建立了一个软件可靠性模型, 然后运用最大似然估计法进行参数估计, 找出了一个比较合理的估计初始错误个数和指数系数的方法。
关键词：软件可靠性,JM模型,参数估计,故障率
参考文献
[1]尚珊珊, 赵轶群.软件可靠性综述[J].软件导刊, 2006 (8) :3-5.
[2]邹丰忠, 李传湘.软件可靠性模型理论分析[J].武汉水利电力大学学报, 1998, 31 (2) :76-80.
[3]朱东鸣.MATHEMATICA软件建立软件可靠性模型[J].成都大学学报 (自然科学版) , 1996, 15 (2) :32-37.
[4]徐仁佐, 谢云, 郑人杰.软件可靠性模型及应用[M].北京:清华大学出版社, 1994.
[5]谭民, 疏松桂.软件可靠性模型的发展[J].计算机学报, 1990 (5) :382-390.
[6]王峰, 陈杰, 喻林.计算机软件测试[M].北京:机械工业出版社, 2004.
非参数估计法 篇6
非平稳信源是现实世界中较为普遍存在的一类信源, 因其统计特性随着时间的推移而不断的变化, 所以成为信源编码领域中较难处理的一类信源。目前, 绝大多数的分布式信源编码系统都是建立在两个相关信源的相关性信息的统计特性是非时变的基础之上, 通过充分利用斯理篇-伍夫[1]定理的理论极限来实现信源的无损或是一定保真度内的有损压缩。然而在许多实际的分布式信源编码的应用中, 编码器无法预先获知相关信源的交叉概率, 因此编码器需要首先估计出相关性信息以便确定编码速率, 最终达到压缩信源的目的。
在分布式信源编码的框架下, 针对相关性统计特性的估值问题, 许多学者提出了不同的实现方案, 总体上分为两大类:离线估计方案和在线估计方案。目前绝大多数是基于非现实的针对平稳信源的离线估计方案, 即假设边信息或者部分边信息在编码器端可见或者是原信源的部分信息在译码器端可见[2,3,4,5]。同时由于斯理篇-伍夫编码与信道编码之间的密切关系, 这些信源的参数离线估计方案也经常被用来进行信道噪声的估计[6]。一种针对像素域与转换域的视频编码相关性噪声建模应用则进行了相关信息统计特性的在线估计[7]。
然而现实世界中更普遍存在的随机过程是非平稳的, 因此深入研究非平稳信源的参数估计问题具有重要的意义。本文探索了一种基于隐马尔科夫模型的非平稳信源参数估计方法。首先建立模型, 将非平稳信源与隐马尔科夫模型对应起来;接着对隐马尔科夫模型的前向算法和后向算法进行改进并利用改进后的算法计算相应的变量;最后根据非平稳信源每一时刻偏移概率的实际意义, 利用算法得到的结果对非平稳信源进行参数估计。
1 基于HMM的非平稳信源参数估计系统框架
1.1 隐马尔科夫模型的构造
HMM是在Markov链的基础上发展起来的, 由于实际问题比Markov链模型所描述的更为复杂, 观察到的事件并不是与状态一一对应的, 而是通过一组概率分布相联系, 这样的模型被称为HMM[8,9,10]。它是一个双重随机过程, 其中一个描述的是状态之间的转移;另一个描述的是状态和观察符号之间的统计对应关系。观察符号与状态之间并没有一一对应的关系, 因此, 只能通过观察符号感知到状态的存在及其特性。
隐马尔科夫模型作为一种统计工具, 其模型本身也是由一系列参数来刻画其特性的。一阶隐马尔科夫模型主要由隐状态数目N、观察符号数目M、隐状态转移概率矩阵A、观察符号概率矩阵B以及隐状态初始概率分布向量π来描述其特性。由于本文主要研究二进制非平稳信源的统计特性, 因此观察符号数目M选择为2;为了更好地逼近非平稳信源的统计特性同时尽量降低计算量, 将隐状态数目N选为27;由于高斯分布更加接近现实生活中的随机过程, 因此状态转移概率矩阵A选为服从高斯分布的概率矩阵;为了降低计算复杂度, 观察符号概率矩阵B选为与隐状态成正比的整数;隐状态初始概率分布向量π选为等概分布。需要指出的是本文的隐马尔科夫模型不同于传统的隐马尔科夫模型, 而是经过了一定程度的变形, 目的是使变形后的模型参数能够适应本文的研究问题。
1.2 非平稳信源参数估计原理
本文首先将非平稳信源的输出对应为隐马尔科夫模型的可观测序列;其次对选取的隐马尔科夫模型进行适当的变形, 主要是调整隐状态转移概率矩阵A、观察符号概率矩阵B以及隐状态初始概率分布向量π以解决模型参数过小所导致的计算精度问题;接着改进前向算法和后向算法, 以确保算法适用于本文的研究问题;最后对观测序列的时变统计特性进行实时估计。其系统框架如图1所示。设o={ot}Tt=1代表二进制的非平稳信源, 即ot∈{0, 1};p={pt}Tt=1代表该信源的偏移概率, 其中pt=Pr{ot=1}, 即pt代表某一时刻观察值ot等于1的局部偏移概率。本文的目的就是根据观察值序列, 结合选定的隐马尔科夫模型, 估计出决定该观察值序列的偏移概率向量p。
2 算法改进
2.1 隐马尔科夫模型的参数调整
为了解决隐马尔科夫模型参数过小所导致的计算精度问题, 本文主要调整了隐状态转移概率矩阵A、观察符号概率矩阵B以及隐状态初始概率分布向量π, 这也是本文构造的隐马尔科夫模型与传统的隐马尔科夫模型的不同之处。本文的隐状态转移概率矩阵A是由一个变形的高斯分布所产生的对称矩阵, 决定该高斯分布的主要参数为其方差σ。其中状态之间的转换概率根据状态之间的汉明距离决定, 即距离越小, 转换概率越大;距离越大, 转换概率越小。观察符号概率矩阵B则是根据当前时刻的输出值, 选取一个对应的与所处状态成正比的整数值。而隐状态初始概率分布向量π的各个分量固定为1, 以确保每一种隐状态在初始时刻具有同等的重要性。
2.2 前向算法
在采用HMM来完成实际的各项课题研究时, 前向算法用来进行输出概率的计算问题。即给定隐马尔科夫模型λ= (π, A, B) , 和观察值序列o={ot}Tt=1, 如何有效地计算观察序列O对隐马尔科夫模型λ的输出概率P (o∣λ) 。然而本文利用前向算法的目的是计算每一时刻每一种状态所对应的前向变量αt (j) , 然后再计算同一时刻所有状态前向变量的统计平均值作为该时刻使用前向算法的估值结果, 记为p Estimate Ft。其中前向变量αt (j) 表示隐马尔科夫模型λ在t时刻, 处于状态sj时观察到部分序列 (o1, o2, …, ot) 的概率。即αt (j) =P (o1, o2, L, ot, zt=sjθ) , j=1, 2, …, N, t=1, 2, …, T。具体的算法实现步骤如下:
(1) 初始化:α1 (i) =πibi, o (1) , 其中i=1, 2, …, N。
2.3 后向算法
给定隐马尔科夫模型λ= (π, A, B) , 和观察值序列o={ot}Tt=1, 针对隐马尔科夫模型的后向算法, 本文提出了两种实现方法, 如下所示:
第一种方案的实现步骤如下:
(1) 初始化:βT (i) =1, 其中i=1, 2, …, N。
第二种方案的实现步骤如下:
(1) 初始化:βT (i) =1, 其中i=1, 2, …, N。
2.4 混合算法
当计算出每一时刻的前向变量αt (j) 和后向变量βt (i) 以后, 可以通过综合利用这两个变量计算混合变量γt (i) , 然后再利用γt (i) 进行参数估计, 记为p Estimate Mt。具体的实现步骤如下:
(1) 计算γt (i) :γt (i) =αt (i) ×βt (i) , 其中i=1, 2, …, N;t=1, 2, …, T。
3 实验结果与讨论
在利用本文提出的方案进行非平稳信源的参数估计时, 将非平稳信源序列长度设为T=212, 图2-图5分别显示了使用本文实现的三种算法对一种非平稳信源进行参数估计的效果图。其中p表示原始的非平稳信源;p Estimate F表示使用前向算法进行参数估计得到的估值结果;p Estimate B表示使用后向算法进行参数估计得到的估值结果;p Estimate M表示使用混合算法进行参数估计得到的估值结果。
分析图2-图5可以得出:对隐马尔科夫模型的前向算法和后向算法进行适当地改进之后可以用来进行非平稳信源统计特性的估计;同时可以看出混合算法得出的估值结果最为逼近原始的非平稳信源, 这一点也可以从三种算法的估值结果与原始信源之间的均方误差中得出, 前向算法得出的均方误差为0.005148, 后向算法和混合算法分别是0.004553和0.002210, 这也说明通过综合利用前向算法和后向算法可以更好地逼近原始信源的变化规律。图6和图7分别显示了使用两种不同实现方案的后向算法对同一种非平稳信源进行参数估计的效果图, 其中的变量说明同如2-图5。
分析图6和图7可知:两种不同的后向算法都可以用来进行非平稳信源的统计特性估计。本文通过大量的实验分析两者之间的优劣, 最终得出两种方案对不同的非平稳信源进行参数估计时的效果差距不大, 各有优点, 都能用于逼近原始信源的变化规律。表1记录了对十五种随机产生的非平稳信源进行参数估计的结果, 其中p表示原始的非平稳信源, Error F列表示p Estimate F与p的均方误差;Error B列表示p Estimate B与p的均方误差;Error M列表示p Estimate M与p的均方误差。
分析表1的数据可以得出与图2-图5相同的结论:即隐马尔科夫模型可以用来估计非平稳信源不断变化的统计特性;另外相比前向算法和后向算法, 混合算法得出的估值结果能更好地逼近原始信源的变化规律, 为后续的数据压缩提供了实现基础。
4 结语
针对非平稳信源不断变化的统计特性, 本文提出了一种基于HMM的参数估计方案。通过构造合适的HMM, 改进前向算法和后向算法, 并在此基础上提出估值效果更优的混合算法, 最后对随机产生的大量非平稳信源进行验证。实验结果表明应用本文提出的参数估计方案可以较准确地掌握非平稳信源的变化规律, 为进一步实现该类信源的压缩提供了实现基础。
参考文献
[1]Slepian D, Wolf J.Noiseless coding of correlated information sources[J].IEEE Trans.Information Theory, 1913, 19:471-480.
[2]Cheung N M, Wang H S, Ortega A.Sampling-based correlation estimation for distributed source coding under rate and complexity constraints[J].IEEE Trans.Image Processing, 2008, 17 (11) :2122-2137.
[3]Fang Y, Jeong J C.Correlation parameter estimation for LDPC based Slepian-Wolf coding[J].IEEE Commun.Lett., 2009, 13 (1) :67-76.
[4]Fang Y.Crossover probability estimation using mean-intrinsic-LLR of LDPC syndrome[J].IEEE Commun.Lett., 2009, 13 (9) :679-681.
[5]Varodayan D.Aaron A, Girod B.Rate-adaptive codes for distributed source coding[J].Signal Process., 2006, 86 (11) :3123-3130.
[6]Cui L, Wang S, Cheng S, et al.Noise adaptive ldpc decoding using particle filters[C]//Proc.CISS 2009, Baltimore, MD, 2009.
[7]Brites C, Pereira F.Correlation noise modeling for efficient pixel and transform domain Wyner-Ziv video coding[J].IEEE Trans.Circuits Syst.Video Technol., 2008, 18 (9) :1177-1190.
[8]Lawrence R Rabiner.A Tutorial on Hidden Markov Models and Selected Applications in Speech Recognition[J].Proceedings of the IEEE, 1989, 77 (2) :257-286.
[9]Wang Dong, Huang Lin.HMM Based Distributed Arithmetic Coding and Its Application in Image Coding[J].ICIMT.2010, 2:397-400.
非参数估计法 篇7
随着计算机视觉技术的发展,基于图像处理的生物量估计成为智能化水产养殖的研究热点[4]。现有的基于图像处理的数量估计方法分为两种:一种是逐个计数[5,6,7,8],使生物群经过摄像区域的通道,通过图像处理实现自动计数,该方法计数准确率高,但需要将生物群进行转移,对生物个体有损伤且操作繁琐;另一种是采样估计[9],在一些局部区域内进行视频录制,将局部生物量作为样本数据,通过概率模型估计生物数量,该方法操作简单、适应性强,但视频采样的使用率不高,采样方法以捕捉和人工观测居多[10,11,12],存在人力成本高、主观性强的问题。
本研究提出一种基于计算机视觉和非参数估计的数量估计方法,首先使用摄像装置在养殖池内分水层采样,再对采集视频进行图像处理,然后使用轮廓提取得出样本数据,最后使用合适的概率模型对池内蟹苗数量进行估计。
1 视频采集方案
1.1 环境与硬件
试验在圆形养殖池(直径1.2 m,水深0.53 m)内进行。拍摄对象为溞状幼体一期、体长1~2 mm的蟹苗。硬件设备包括转动装置和摄像装置,转动装置转速0.1 r/s,转动杆半径0.6 m。摄像头采用Basler acA1600-20gm工业摄像头,其下端蟹苗通道的入口高度为0.06 m,拍摄区域高度为0.015 m,拍摄范围0.1 m×0.1 m,通道底部为相同尺寸的红外背光源。
1.2 视频采集
视频采集装置结构如图1所示,转动电机固定于养殖池上方,摄像装置上端固定在转动装置的转动杆上,距离电机0.5 m处,下端为透明的蟹苗通道,摄像机放置在蟹苗通道上方的密封空间里。转动装置带动摄像装置在养殖池内匀速转动,蟹苗通过蟹苗通道进入拍摄区后流出通道。视频采集过程将养殖池划分为5个水层,以地面为起始点向上延伸,设置采集高度分为0.08 m、0.18 m、0.28 m、0.38 m和0.48 m。摄像装置通过调整固定点的位置来改变视频采集的高度,以不同高度条件下采集的视频作为对应水层的样本视频。
1.3 图像预处理
采用基于混合Gauss模型的背景建模方法(GMM)[13]获取背景图像,然后将当前图像与背景图像相减得到前景图像。图2为原图与前景提取图。
GMM算法对图像背景的变化有较好的适应性,可以排除部分水下环境杂质的影响。为进一步去除图像中的杂质,本文对所得的前景图像提取轮廓,将轮廓周长小于正常蟹苗图像周长的像素块去除,处理结果二值化后如图2b,其中白色连通域表示蟹苗。
2 蟹苗数量估计
蟹苗数量估计流程如图3所示。首先将采集的图像序列进行预处理得到二值图,再以轮廓提取方法统计图像序列中的蟹苗总数,然后将各序列的蟹苗总数作为对应水层的样本数据,通过核密度估计,得出概率密度函数,最后结合概率密度函数和样本数据,估算出养殖池内的蟹苗数量。
2.1 轮廓提取
通过图像预处理后得到前景图像,然后使用轮廓提取即可得到图像中的蟹苗个数。但图像中常有蟹苗粘连现象,导致提取的轮廓数低于实际值。因此,先对原图使用分水岭算法,得到蟹苗对应的连通域,然后使用形态学滤波的方法,对连通域进行腐蚀运算,将粘连的连通域分离。由图4中的矩形圈定的区域可以看出,在进行分水岭算法处理之后,部分粘连蟹苗成功分离;经腐蚀处理后,所有粘连蟹苗都可以成功分离。对所有粘连蟹苗分离之后,再对结果图提取所有轮廓,即可得到图像中的蟹苗数。
2.2 核密度估计
核密度估计是利用某一种核函数的线性组合来估计概率密度分布的方法,也就是通过一定区间范围内的各样本点密度的均值来估计总体密度函数[15,16,17,18]。该方法通过样本数据集可以产生一个渐近无偏的密度估计,设定一组样本S={xi}i=1…N的概率密度函数为p(x),以x为中心点作一个超立方体。超立方体的边长为h,体积hD,包含样本数为n。
定义窗函数如式(1):
超立方体中的样本数如式(2):
则概率密度估计如式(3):
2.3 数量估计
通过对视频数据的处理得出每个水层的样本数据。将这些样本数据转化为对应水层的概率值,先通过核函数估计,将单个水层的概率值以核函数的分布规律分散于以样本为中心的整个空间,所得的分布曲线则为概率密度曲线,然后将多个水层对应的多个曲线叠加,得到整个空间的概率密度曲线。其中,核函数选取高斯函数,带宽为1。图5为实验中得出的概率密度函数图,该图直观地显示出养殖池内蟹苗的分布情况。
对概率密度函数积分求概率值时,积分区间半径的选择非常重要。区间半径过大或过小都会影响到估计值的准确性。估计值变化曲线如图6所示。
因此,需要不断调整区间半径,使得估计值与实际值的差值减小。图6显示的是在不同区间半径下,3组实验的概率密度函数得出的估计值变化曲线与实际值直线(y=96 000)的相交情况。结果显示,3组实验的实验对象和实验环境相同,实际蟹苗数量约为96 000只。
估计值曲线与实际值直线(y=96 000)相交时,对应的区间半径最佳。因此,可以通过统计不同密度条件下估计值曲线与实际值直线的交点,得到最佳区间半径的大致范围。表1为不同密度条件下,所得的最佳区间半径。
由表1可以得到,最佳区间半径都在2~3之间,因此,本文选取区间半径为2和3估计值的平均值,作为池内蟹苗的估计量。其中,估计量的计算公式如式(4):
式中:f(x)为核函数估计所得的概率密度函数;N(x)为水层高度对应的蟹苗数量统计值;r1和r2为区间半径。对f(x)在区间内积分即可得出对应的概率值,将该水层的统计值除以概率值,即得出蟹苗总量的估计值。
3 结果
在视频录制之后,将蟹苗分装在多个桶中,然后对每一桶使用1 L的烧杯进行抽样统计,统计得出的密度均值作为实际蟹苗密度。实验分别在蟹苗平均密度为120、140和160只/L的养殖池内进行视频采集,通过图像预处理和轮廓提取后,得出视频中蟹苗的个数。不同采样高度下视频中的蟹苗总数如表2所示。
对表2数据通过核密度估计得出概率密度函数后,结合表2数据和本文统计得出的最佳区间半径范围,即可估算出池内蟹苗数量。表3为不同密度条件下,实际值与估计值的对比。其中,标准值由实际蟹苗密度乘以总水量而得出。由表3可看出,低密度条件下的估计正确率较高,高密度条件下的估计正确率较低。
在不同蟹苗密度条件下分别获取的100帧图像,使用背景建模和轮廓提取的方法进行图像处理,得出累计的蟹苗个数作为估计个数。图像处理后,在图像中标记出目标个体,然后通过人工判别得出实际个数。表4中显示的是在不同密度条件下,估计个数和实际个数的对比。
4 讨论
4.1 不同密度条件下正确率的对比
低密度条件下的估计正确率较高。其原因是,在该条件下图像中蟹苗个体间隔较大,导致图像处理效果较好。而在高密度条件下,蟹苗个体之间会发生遮盖,导致图像处理结果较差。并且,高密度条件下养殖池内的杂质增多,对图像处理影响较大,导致估计正确率较低。图像中判别蟹苗的正确率直接影响了蟹苗总量估计的正确率。因此,针对高密度条件下的蟹苗个体识别需要进一步的改善。
4.2 计数方法的对比
王硕等[5]针对大菱鲆鱼苗的计数,提出了一种基于曲线演化的方法来分离粘合图像,从而得出图像中的鱼苗个数;朱从容等[8]利用计算前景像素块和鱼苗个数的关系对鱼苗进行计数;王文静等[7]针对多帧图像提出一种连续实时在线计数方法,使鱼苗匀速流过拍摄区域并进行鱼苗计数。上述方法研究对象都需要被转移至环境简单的图像采集区域,本方法主要针对实际育苗环境下的蟹苗计数,使用背景建模和杂质去除来消除环境的干扰,这是一种从实验环境到实际环境的探索。并且本文的计数方法为轮廓提取后计数,不受对象的大小限制,可适用于混合养殖环境下的计数。
在实验中,本方法在低密度条件下处理良好,但在高密度条件下处理效果不佳,图像中常出现絮状物与蟹苗混杂的情况,需要对两者分割处理后才能得到更准确的计数。因此,本方法在排除环境杂质的图像处理上需要进一步改善。
4.3 数量估计方法的对比
常用的数量估计方法有均匀分布估计法[19]和直方图估计法[20]。前者使用比较普遍,它以多个局部密度的均值代替整体密度来获取估计值;后者则是通过对各个区域的生物分布情况进行统计,得出整体概率密度分布,从而获得估计值。本文对同组数据分别使用这两种方法和本文的方法进行数量估计,其中,均匀分布估计法正确率7342%,直方图估计法正确率80.57%,本文估计方法正确率82.14%。可以看出,均匀分布估计法虽然实用性很强,但在分布不均匀的情况下,正确率不高;而直方图估计法,由于其概率密度曲线不平滑,对未采样的区域估算误差较大,使得估计正确率较低。
4.4 本估计方法的推广
本文方法的研究对象除了蟹苗,也可以是鱼苗、虾苗。其中,鱼苗个体比蟹苗的大,且不容易与杂质混合,图像比较好处理;虾苗的活动习惯和蟹苗相近,其分布较为均匀,适合使用概率模型统计。因此,本文的方法也可推广到鱼苗和虾苗的数量估计。
4结论
提出了一种基于计算机视觉和非参数估计的蟹苗数量估算方法。该方法可以比较准确地估计蟹苗数量。但在高密度的情况下,蟹苗之间有遮盖现象,且池内杂质增多,会产生图像处理结果不佳的结果,因此,该图像处理方法有待进一步提高。受实验环境的限制,本方法估测范围偏小,在以后的研究中,将针对大型养殖池的蟹苗数量估计进行更多的探索。
摘要：养殖池内蟹苗的数量估计在蟹苗养殖中有着重要意义。但现有的数量估计方法操作复杂且实用性不强,因此提出一种基于计算机视觉和非参数估计的蟹苗数量估算方法。首先在养殖池内分水层采集视频,通过背景建模得到前景图像,并使用分水岭算法和轮廓提取得到视频中的蟹苗数量,并作为样本数据,然后通过核密度估计得出概率密度函数,最后结合该函数和样本数据估算出池内蟹苗数量。结果表明,该方法对于容积约为1 000 L、蟹苗密度100~160只/L的小型蟹苗养殖池,估算蟹苗数量的平均正确率为82.14%。研究表明,采用该方法不仅可以解决采集视频过程的操作繁琐、幼苗转移的问题,而且能够避免图像处理过程中部分背景杂质的干扰。该方法还可以推广到虾苗和鱼苗等生物的幼苗估计,具有良好的通用性和可行性。
非参数估计法 篇8
由于电子装置、半导体器件等非线性负荷的广泛使用,电力系统谐波污染日益严重。为了有效地治理谐波污染,就必须对电网谐波参数进行准确检测和分析。快速傅里叶变换(FFT)是电网谐波参数估计最直接、最有效的方法。但非同步采样时存在的频谱泄露和栅栏效应,会导致谐波参数估计的精确度下降[1]。
针对FFT存在的问题,国内外学者提出一系列的高精度估计方法。文献[2-5]采用插值FFT算法在一定程度上提高了谐波参数计算的精度,但是随着插值拟合函数的阶次增高及谐波次数的增多,谐波参数估计的计算量会大大增加。文献[6-7]中提出的能量重心法是另一种通用的参数估计方法,其原理是根据各种窗函数离散频率的能量重心无穷逼近坐标原点或在原点附近的特点进行参数估计。本文通过理论推导和数值试验,讨论了影响能量重心法精度的因素,并且分析了能量重心法的抗噪性能。
1 算法原理
设一个频率为f0、幅值为A、初相为单一频率信号x(t),在经过采样率为fs的模数转换后得到离散信号
对信号x(n)加窗截断,则加窗信号x(n)w(n)的频谱为
其中:W(k)是窗函数w(n)的频谱;。
设m为幅值最大谱线序列号,d为由非同步采样引起的频率偏移量,当最高谱线在次高谱线左侧时,k0=m+。只考虑在正频率f0附近的频谱有
考虑到窗函数的对称性,式(3)右边第一项为零,化简整理后有
如果选择合适的窗函数,使得绝大部分能量集中在主瓣内,则用主瓣能量代替总能量带来的误差就很小;另一方面,加窗函数只能减小“长泄漏”,即“旁瓣干扰”的影响,分析多频信号的电力谐波,考虑“长泄漏”的影响是困难的,因此,也只能用窗函数的主瓣能量来替代其总能量。设主瓣内谱线数目为B2k,则式(3)可化简为
则信号的频率和相位的估计公式分别为
其中,φm为信号幅值最大谱线的相位。
另一方面,由Parseval定理[8],并考虑到加窗后由于窗函数的加权作用有
其中,Ke为窗函数的能量恢复系数[9]。
考虑到能量主要集中在主瓣内,式(8)可近似化简为
也就是
式(6)、式(7)、式(10)即是最高谱线在次高谱线左侧时能量重心法的频率、相位和幅值的估计公式。同理也可以推出最高谱线在次高谱线右侧时的估计公式。
能量重心法的实质是通过主瓣内的谱线来近似计算信号的能量,进而得到各参量的估计公式。其估计精度取决于所加窗函数在主瓣内能量集中度。
对于单频信号,参与校正的谱线数越多,则能量重心的估计精度也越高,故采用能量重心法估计单频信号参数时,校正点数甚至可以取到旁瓣。但电力系统中的谐波是多频率成分信号,由于谱间泄漏将造成各次谐波之间的谱间干扰,理论上FFT各谱线的值是各信号分量单独经FFT后的矢量和,因此,用增加参与校正的谱线数来提高估计精度的方法是不可行的,只能采用能量主要集中在主瓣的窗函数来加权截断被分析信号。这样,一方面可减小“长泄漏”的影响,即使FFT的谱线主要是一个频率信号的量值;另一方面可提高用主瓣能量代替总能量的精度。进而提高了用kB条谱线分析各谐波参数的分析精度。故对于多频信号,应只取窗函数主瓣内的谱线kB进行参数估计。
2 仿真实验
本文采用的仿真软件为Matlab R2010a。仿真采用的信号模型为
式中:基波频率f0为50.4 Hz;FFT时,一般要求采样点为2N个,本文取210=1 024个点;IEC标准规定采样时间为10个周波(0.2 s)[10];故采样频率为1024/0.2=5 120 Hz;仿真的电力信号的相关参数如表1所示。
2.1 多频信号参数估计
采用Blackman窗截断。考虑到Blackman窗的主瓣内有6条谱线,故参与估计的谱线为6;Ke为3.779 3。利用能量重心法,高精度估计出基波频率、幅值和相位为
图1、图2为文献[11]的双峰谱线插值法和本文的能量重心法的幅值和相位估计的误差曲线。可见两种方法都有很高的分析精度,且能量重心无论是幅值还是相位都具有更好的分析结果。用同一台计算机和同一个仿真软件进行计算,计算耗时上,双峰谱线插值法为0.015 830 s,能量重心法为0.006 488 s,可见能量重心法所需时间仅为插值法的一半。
2.2 窗函数的影响
不同的窗函数在主瓣内的能量集中度是不一样的。当采用不同的窗函数时,能量重心法的估计结果的精度也是不一样的。图3、图4是分别采用Hanning窗、Kaiser窗(β=7.856)以及Blackman窗时,能量重心法幅值和相位估计的误差曲线。
分析可知,Blackman窗对于基波参数的估计精度比Hanning窗和Kaiser窗都高出一个数量级,这是因为这三种窗函数中,Blackman窗的主瓣内功率相对于总功率比重最高,也即能量集中度最高。可见仿真结果是和理论分析一致的:窗函数主瓣能量集中程度越高,则参数估计精度越高。
2.3 高斯白噪声对估计精度的影响
噪声的存在会使频谱的幅值和相位值发生变化,影响信号谐波参数估计的准确性[12]。结合式(11)的信号模型,对高斯白噪声存在情况下的谐波参数估计进行仿真实验。
仍采用Blackman窗,参数估计的谱线取6。
图5、图6是加入基波信噪比为50 d B的噪声时,幅值和相位估计的误差曲线。从曲线可知,噪声的存在降低了能量重心的精度,尤其是对相位的估计精度。对于含量比较小的16次、20次以及21次谐波,三种分量对应的信噪分别为-21.29 d B、-42.85 d B和-36.83 d B,这些弱信号分量几乎被噪声淹没了,参数的估计误差很大。
图7、图8是不同信噪比下,能量重心法和双谱线插值法对基波参数估计的误差曲线。
由图7、图8知,无论是信噪比高低,双峰谱线插值法和能量重心法对基波参数的估计精度都很高。对于幅值的估计,能量重心法是优于插值法的:当信噪比增大时,插值法的幅值估计误差基本恒定[13],但是能量重心法幅值的估计误差迅速下降。对于相位估计,能量重心法精度较双峰谱线插值法提高了1~2个数量级。可见能量重心法具有良好的抗噪性能[14]。
3 结论
(1)能量重心法可以大大提高电力谐波参数估计的精度,算法简单,很容易在数字系统上实现,计算速度快,有利于提高系统的实时性。
(2)对于电力系统谐波信号,由于存在频谱泄露现象,考虑到信号的能量主要集中在主瓣内,故仅用主瓣频谱来对就可以实现对参数的高精度估计。
非参数估计法 篇9
电力系统动态仿真中所使用的元件模型,特别是电力系统动态元件——发电机参数的准确性对精确计算和分析电力系统的行为具有重要的意义。而现场往往缺少实际参数,在计算中只能查用工厂或手册的典型数据。由于这些数据未计及涡流、磁滞、饱和等实际运行工况的影响,使得计算结果与实际情况不符,严重影响了仿真计算的准确度和可信度。因此,如何准确估计发电机的参数就成为电力系统动态仿真及稳定分析的重要内容之一。
关于发电机参数的辨识已经有很多方法,概括起来可以大致分为2类:一类是频域响应法[1,2],该方法在计算方法上比较成熟,算法稳定性好,但对输入扰动信号的波形、幅值大小及其相关性要求严格,存在不能够反映动态过程中参数非线性变化等问题;另一类是时域辨识法[3,4,5,6],其最大优点是能够将饱和、涡流、旋转等因素自然地包含在参数估计值中,缺点是存在参数辨识结果不稳定的问题。时域辨识方法主要有最小二乘法和卡尔曼滤波法,这2种方法一方面存在局部最小或收敛性的问题,另一方面采用模型线性化方法会产生固有误差[7]。
近年来,相量测量单元(PMU)技术的快速发展和广域测量系统的推广应用为发电机参数在线辨识与估计提供了新手段。PMU量测数据直接反映物理系统的响应,计及了发电机实际运行工况,一旦辨识成功,将具有较高的可信度;而且所采用信号来自物理系统扰动,不需要人为进行试验,大大节省了人力、物力。基于功角测量的模型及参数辨识方法既可以采用时域辨识法[8,9],也可以采用在线频率响应方法[10]。实际应用中由于测量设备本身的缺陷及数据传输过程中随机扰动的影响,使得测量结果中不可避免地存在误差和坏数据,这会影响参数辨识结果的准确性,因此,需要研究效率高、抗差能力强的估计模型和算法。本文将一种具有很强抗差能力的最大指数平方估计模型[11]引入到参数估计中,并采用轨迹灵敏度进行求解。本文方法可应用于在线的模型校正以提高仿真计算结果的可信度。
1 轨迹灵敏度
轨迹灵敏度法描述了由于基础参数或者初始条件的改变导致的系统运动轨迹变化的线性化关系。轨迹灵敏度方法已经被成功应用于电力系统安全稳定评估、电压紧急控制[12,13]、预测状态变量的轨迹变化[14]、模型预测控制[14,15]、电力系统动态安全控制[5]等方面。本文把轨迹灵敏度方法引入到参数辨识中,其基本原理介绍如下:
假设电力系统动态过程中状态变量x与代数变量y可以分别表示为:
式中:λ为模型参数,本文主要考虑发电机的电抗参数xd,xq,xd′,xd″,xq′,xq″和时间常数Td0′,Tq0′,Tq0″,Td0″,Tj;x0为初始条件。
对式(1)模型参数进行一阶Taylor级数展开,得到:
对于y有同样的形式:
式(2)和式(3)中的xλ和yλ分别为x和y 关于参数λ的轨迹灵敏度。
设系统动态方程如下:
求解动态方程获得的动态响应过程对各参数轨迹灵敏度的具体表达式如下:
式中:xλ的初始条件为:xλ(0)=0;yλ的初始条件可由下式获得:0=gy(0)yλ(0)+gλ(0)。
对式(5)、式(6)求解即可获得输出结果对参数的一阶轨迹灵敏度[16]。
同理,对式(1)进行二阶Taylor级数展开时,可以获得输出结果对参数的二阶轨迹灵敏度。
2 估计模型
发电机参数较多,全部估计无法保证可观测性。实际上,在这些参数中只有部分对动态稳定性影响较大,称之为主导参数。因此,首先利用轨迹灵敏度寻找出对仿真结果产生较大影响的主导参数集作为待估计参数集,然后将动态过程中的PMU量测值与模型仿真的结果进行比较,根据两者的不同对主导参数进行调整,以使两者获得最佳匹配,即得到参数的估计值。
由于量测设备本身的缺陷及数据传输过程中随机扰动的影响,使得量测结果中不可避免地存在坏数据,这使得参数辨识过程受到极大影响。因此,本文采用了一种具有抗差能力的估计模型[11]:
式中:ym(ti)和y(ti)(i=1,2,…,N)分别为PMU实测数据及拟合数据,对应于发电机的机端电压及功角。
满足一阶轨迹灵敏度对应的等式约束为:
y(t)=φ(x0,t,λ0)+yλ(t,λ0)Δλ (8)
目标函数式(7)具有很好的性质,当存在坏数据时,坏数据对目标函数的贡献接近于0,因此坏数据对估计结果的影响被自动压缩。同时可以证明,当所有量测的误差都较小时,该估计模型等价于最小二乘法。
将等式约束式(8)代入目标函数式(7)中,可以将有约束优化问题转化为无约束优化问题。令k=1,然后采用牛顿法根据Taylor级数将方程线性化,迭代搜索求得最优的Δλ,对原参数λ更新得到λ=λ+Δλ,然后由式(5)、式(6)计算新的轨迹灵敏度,并令k=k+1,进行再次优化。循环迭代的终止条件可以设为k≥Nmax或前后2次参数的改变量max Δλi<ε[17],ε为一极小数。
另外,还可以采用二阶轨迹灵敏度进行估计,只需把一阶轨迹灵敏度对应的等式约束式(8)改成二阶轨迹灵敏度对应的等式约束式(9)即可。
由于Taylor级数展开式的阶数越高,拟合误差越小,因此,二阶轨迹灵敏度拟合比一阶拟合具有更高的精度。理论上,要达到同样的精度,二阶拟合的迭代次数要小于一阶拟合的迭代次数。
3 算例分析
以新英格兰4机双区域系统为研究对象,对同步发电机的参数进行辨识,系统结构如图1所示。所有发电机均采用考虑饱和的五阶发电机模型。负荷采用感应电动机并联恒阻抗模型。
故障设置:1 s时,在B6与B7之间的一回线路上,靠近B6的2%处发生三相接地短路故障,1.1 s切除这一回故障线路。
以发电机G1的参数为辨识对象,仿真得到发电机机端电压、电流及励磁电压、转速、功角的变化数据,利用这些数据作为发电机的实测数据,对同步发电机模型进行辨识。
首先分析发电机G1机端电压及发电机G3的功角对G1电抗、时间常数等参数的轨迹灵敏度。G1机端电压对其各电抗参数的轨迹灵敏度见图2。
从图2中可以看出,xd,xq,xd′,xd″轨迹灵敏度大,xq′,xq″灵敏度小,且具有明显的时域特征,故障刚发生时刻有一定的作用,但是这种作用随时间的增加而迅速衰减,这与发电机的暂态、次暂态和稳态对应。
G1机端电压对其各时间常数的轨迹灵敏度如图3所示。可以看出,Td0′,Tq0″,Tj作用最明显,其他参数的作用不明显。它们的作用与电抗类似,也与发电机的暂态、次暂态和稳态对应。响应为G3功角时确定出的主导参数与响应为G1机端电压时的结果一致,限于篇幅,这里不再给出详细结果。
通过上面的分析可以看到,xd,xq,xd′,xd″,Td0′,Tq0″,Tj这7个参数的轨迹灵敏度较大,对响应的影响较大,因此,下面将以这7个参数作为需要辨识的参数。
假设在1 s时, B5与B6之间的一回线路上,靠近B5的2%处发生三相接地短路故障,1.1 s切除这一回故障线路,以仿真得到的发电机机端电压、功角的变化数据作为实测数据对这7个参数进行辨识。
对这7个参数分别叠加上幅度为20%的偏差,辨识结果如表1所示。可见,假如量测没有误差,本文方法得到的估计值精度很高。
对量测真值数据加上标准偏差为0.03的正态分布误差,量测真值及叠加扰动后的数据如图4所示。
利用叠加扰动后的数据作为发电机的实测数据对同步发电机模型进行辨识,结果如表2和表3所示。从表2和表3可以看出,估计误差较小,表明最大指数平方法估计参数鲁棒性好。
本文还采用二阶轨迹灵敏度对发电机参数进行了估计。表4给出了在参数偏差为40%时二阶轨迹灵敏度的估计结果。
对比表4和表3的结果可以看出, 基于二阶轨迹灵敏度的估计结果与基于一阶轨迹灵敏度的估计结果基本一样。但是一阶估计需要8次迭代,而二阶估计大约需要5次迭代,这是因为二阶灵敏度的拟合误差较小。
进一步加大白噪声的标准偏差到0.08,辨识结果如附录A表A1所示。可见,在没有坏数据的情况下,最大指数平方法和最小二乘法的估计结果精度基本一致。
随着噪声的进一步增大,2种方法的辨识结果均开始变差,本文方法略优于最小二乘法。这是因为随着噪声的增加,指数型目标函数会自动压缩突变相对较大的噪声数据的影响。噪声标准差加大到0.15时的辨识结果如附录A表A2所示。
由于PMU量测可能存在坏数据,这会影响参数估计结果。附录A表A3给出了量测数据中存在坏数据时的辨识结果。附录A表A3的结果表明,最大指数平方法可以自动剔除坏数据的影响,而最小二乘法的估计结果明显恶化,需要进行坏数据辨识后再进行估计。进一步的实验表明,将坏数据的比例增加到5%时,最大指数平方法仍具有很好的抗差能力,而最小二乘法的估计结果迅速恶化。
对以上数据的分析可以看出,采用轨迹灵敏度的估计方法是合理有效的。另外,由于采用了最大指数平方估计模型,它克服了对坏数据的敏感性问题,而且在参数变化的一定范围内估计结果的精度还是很好的。
4 结语
本文提出了一种求解发电机参数的方法——轨迹灵敏度法,同时引入了一种新的参数估计模型——最大指数平方模型,可以自动剔除坏数据的影响,并将其应用于发电机的参数估计。在新英格兰系统上的仿真结果证明了本文方法的有效性,具有应用于在线参数估计的潜力。
附录见本刊网络版(http://aeps.sgepri.sgcc.com.cn/aeps/ch/index.aspx)。