状态空间模型中实际参数估计

  • 状态扩增法
  • 线性状态空间模型的参数估计
  • 利用高斯滤波与平滑的参数估计(非线性模型)
  • 基于粒子滤波与平滑的参数估计
  • 参数的 Rao-Blackwell 化

(参数估计所有内容)

状态扩增法

将参数扩充为状态的一部分,例如有一个含有未知参数 θ\bm{\theta}θ 的非线性模型为:
上述模型可以表述为:
两个参数的动态模型是一致的,重新定义状态为:x~=(xk,θk)\bm{\tilde{x}=(x_k,\theta_k)}x~=(xk​,θk​),则状态空间模型变为:

此时,模型中不含任何未知参数。状态扩增法在整个系统为线性时,可取得良好的效果。

线性状态空间模型的参数估计

假设一个含有未知参数 θ\bm{\theta}θ 的线性高斯状态空间模型为:
其中,
线性高斯模型的能量函数:
能量函数的递归表达式:
式中,vk(θ),Sk(θ)\bm{v_k(\theta),S_k(\theta)}vk​(θ),Sk​(θ) 由只含 θ\bm{\theta}θ 的卡尔曼滤波器给出。
卡尔曼滤波预测过程:
卡尔曼滤波更新过程:

因此,若已知 θ\bm{\theta}θ,先计算 k=0\bm{k=0}k=0 时刻的值 φT=−logp(θ)\bm{\varphi_T=-logp(\theta)}φT​=−logp(θ),在依次执行上述算法直到 k=T\bm{k=T}k=T,可得全局能量函数 φT(θ)\bm{\varphi_T(\theta)}φT​(θ)。
已知能量函数后,可用基于Metropolis-Hastings的MCMC采样,产生后验分布的蒙特卡洛逼近等,得到参数等MAP或ML估计。

计算线性高斯模型中的 Q\bm{\mathcal{Q}}Q:
如上诉所示的线性高斯模型,Q\bm{\mathcal{Q}}Q 的表达式为:
式中,所设计的参数由以 θ(n)\bm{\theta^{(n)}}θ(n) 为参数的 RTS\bm{RTS}RTS 平滑器 给出:

EM 算法对线性空间模型有效。

通过分别求:
上式中各个元素的偏导数 ∂Q(θ,θ(n))/∂θ\bm{\partial\mathcal{Q}(\theta,\theta^{(n)})/\partial\theta}∂Q(θ,θ(n))/∂θ,并分别令其为0,可得如下定理:

线性参数模型参数 Q\bm{\mathcal{Q}}Q 的极大化:
当参数选自模型参数:
时, θ∗\bm{\theta^{*}}θ∗ 的极大值为:
可由下式计算得到:

最终,得到关于初始均值 θ=m0\bm{\theta=m_0}θ=m0​ 的极大值为:

线性状态空间模型的EM算法
令 θ\bm{\theta}θ 包含模型参数 {A,H,Q,R,P0,M0}\bm\{{A,H,Q,R,P_0,M_0\}}{A,H,Q,R,P0​,M0​} 的部分子集,利用下述步骤,可获得参数的极大似然估计值。
首先,假设初始参数为值 θ(0)\bm{\theta^{(0)}}θ(0)。
对于 n=0,1,2,⋅⋅⋅\bm{n=0,1,2,\cdot\cdot\cdot}n=0,1,2,⋅⋅⋅ 执行如下步骤:
1). E\bm{E}E 步骤:
利用 θ(n)\bm{\theta^{(n)}}θ(n) 中当前的参数值进行 RTS\bm{RTS}RTS 平滑器解算,并结合平滑结果计算下式:

也就是上上一个定理线性高斯模型中的Q中的式子。
2). M\bm{M}M 步骤:
结合上一个定理:线性模型参数 Q\bm{\mathcal{Q}}Q 的极大化中的式子,计算新的参数,并将其储存在 θ(n+1)\bm{\theta^{(n+1)}}θ(n+1) 中。

Q(θ,θ(n))\bm{\mathcal{Q(\theta,\theta^{(n)})}}Q(θ,θ(n)) 的表达式也提供了一种简单的梯度计算方法,利用 Fisher\bm{Fisher}Fisher 恒定式计算能量函数的梯度。
Fisher\bm{Fisher}Fisher 恒定式:

利用高斯滤波与平滑的参数估计(非线性模型)

在非线性模型中,可以用相应的非线性卡尔曼滤波器及 RTS\bm{RTS}RTS 平滑器进行参数估计。
考虑下述模型的参数估计问题:
式中,
此时,能量函数可由下述基于高斯滤波的算法得到。

基于高斯滤波的能量函数:
能量函数逼近的迭代表达式为:
式中,vk(θ),Sk(θ)\bm{v_k(\theta),S_k(\theta)}vk​(θ),Sk​(θ) 由只含 θ\bm{\theta}θ 的高斯滤波器给出。
高斯滤波预测过程:

高斯滤波更新过程:

上述的能量函数可直接应用于MCMC采样
当然,还有其他如计算ML或MAP估计中的无梯度最优化。

为计算 EM\bm{EM}EM 算法中所需要的期望值,可用高斯假设密度逼近法(如时矩匹配)对下式进行积分逼近:
式中,

积分逼近得到的 Q\bm{\mathcal{Q}}Q 表达式如下:
即:利用高斯平滑计算 Q\bm{\mathcal{Q}}Q:
对本节中的非线性状态空间模型, Q\bm{\mathcal{Q}}Q 表达式如下:
式中的期望值由高斯平滑器得到。
实际中,高斯平滑器和高斯积分可由: EKF/ERTSS\bm{EKF/ERTSS}EKF/ERTSS,或 sigma\bm{sigma}sigma 点方法 进行逼近。sigma\bm{sigma}sigma 点方法包括:Gauss-Hermite或球面容积积分,以及无迹变换等计算。
也可用上一节中线性模型极大化的方法:对于噪声协方差 Q\bm{\mathcal{Q}}Q 的极值表示为:

对于其他参数的 M\bm{M}M 步骤实现,需依据实际的 f,h\bm{f,h}f,h 而定。若是线性的,可用上述类似方法。
若为非线性, Fisher\bm{Fisher}Fisher 恒等式给出一种计算能量函数的方法。

基于粒子滤波与平滑的参数估计

粒子滤波也可用于逼近参数估计中的边缘似然函数及能量函数的计算,在粒子滤波中,考虑的一般模型为:

式中,θ∈Rd\bm{\theta \in \mathbb{R}^d}θ∈Rd 为未知参数。

算法:基于能量函数逼近的 SIR\bm{SIR}SIR:
参数边缘似然的逼近可以在执行**序列重要性采样(SIR)**算法(粒子滤波)中进行计算,
具体步骤为:
先从先验信息中抽取 N\bm{N}N 个采样点 x0(i)\bm{x_0^{(i)}}x0(i)​,即有:
并令 w0(i)=1/N,andi=1,⋅⋅⋅,N\bm{w_0^{(i)}=1/N}, and i=1,\cdot\cdot\cdot,Nw0(i)​=1/N,andi=1,⋅⋅⋅,N。
对于每个 k=1,⋅⋅⋅,T\bm{k=1,\cdot\cdot\cdot,T}k=1,⋅⋅⋅,T,执行如下操作:
1). 从重要性分布中抽取采样 xk(i)\bm{x_k^{(i)}}xk(i)​,即有:
2). 计算权值:
并计算 p(yk∣y1:k−1,θ)\bm{p(y_k|y_{1:k-1},\theta)}p(yk​∣y1:k−1​,θ) 的估计值:
3). 计算归一化权值:
4). 若有效粒子数过少,则执行重采样

参数边缘似然的逼近值为:
相应的能量函数逼近值为:

能量函数逼近可在执行 SIR\bm{SIR}SIR 算法时进行递归解算。

EM\bm{EM}EM 算法可由粒子平滑实现。Q\bm{\mathcal{Q}}Q 逼近的实际表达式根据具体选用的粒子平滑器而定。若选用后向平滑器,可得以下算法:

算法:利用后向平滑器计算 Q\bm{\mathcal{Q}}Q:
设已知模拟轨迹:
并利用后向平滑器,参数为 θ(n)\bm{\theta^{(n)}}θ(n),并利用下式
可得上述 Q\bm{\mathcal{Q}}Q 式为:

若选用重新加权粒子平滑器,可得如下算法:

算法: 利用重新(调整)加权粒子平滑器计算 Q\bm{\mathcal{Q}}Q;
设粒子集为 {xk(i),k=0,⋅⋅⋅T,i=1,⋅⋅⋅N}\bm{\{x_k^{(i)},k=0,\cdot\cdot\cdot T, i=1,\cdot\cdot\cdot N\}}{xk(i)​,k=0,⋅⋅⋅T,i=1,⋅⋅⋅N},并利用重新加权粒子平滑器计算权值 {wk∣T(i),k=0,⋅⋅⋅T,i=1,⋅⋅⋅N}\bm{\{w_{k|T}^{(i)},k=0,\cdot\cdot\cdot T, i=1,\cdot\cdot\cdot N\}}{wk∣T(i)​,k=0,⋅⋅⋅T,i=1,⋅⋅⋅N},并利用下式
可得上述 Q\bm{\mathcal{Q}}Q 式为:
可利用 Fisher\bm{Fisher}Fisher 特性计算能量函数梯度的逼近值。

参数的 Rao-Blackwell 化

利用Rao-Blackwell化方法可以边缘化状态空间模型中的静态参数。有下面的普通摆模型介绍:
设有:
式中,
量测噪声 R\bm{R}R 的方差未知,只给出 R\bm{R}R 的逆 chi\bm{chi}chi 平方先验分布。
获取 R\bm{R}R 分布参数步骤为:
1). 设此前产生的粒子群为 {wk(i),w0:k(i),i=1,⋅⋅⋅,N}\bm{\{w_k^{(i)},w_{0:k}^{(i)},i=1,\cdot\cdot\cdot,N\}}{wk(i)​,w0:k(i)​,i=1,⋅⋅⋅,N} 该粒子群逼近的状态完全分布为:
结合给定的量测值以及采样状态信息时,R\bm{R}R 的条件分布为:
由此,关于状态及参数的完全分布为:
2). 对重要性分布进行采样:
3). 结合 xk(i)\bm{x_k^{(i)}}xk(i)​ 和此前的量测量,可得当前量测量的似然为:
式中,student\bm{student}student 的 t\bm{t}t 分布参数为:
由此,可以计算出 SIR\bm{SIR}SIR 算法中,下一步的重要性权值,进而可计算有效粒子数,并判断是否需要进行重采样:
4). 结合量测量 y1:k\bm{y_{1:k}}y1:k​ 及状态量 x0:k(i)\bm{x_{0:k}^{(i)}}x0:k(i)​,可进一步计算出 R\bm{R}R 的条件分布:
5). 令 k←k+1\bm{k\leftarrow k+1}k←k+1,跳至步骤 1)。

上述的方法就是 Rao-Blackwell粒子滤波器,其中静态参数 R\bm{R}R 已被边缘化,并通过对 vk(i),RK(i)\bm{v_k^{(i)},R_K^{(i)}}vk(i)​,RK(i)​ 充分统计,结合过去粒子 x0:k(i)\bm{x_{0:k}^{(i)}}x0:k(i)​ 和量测量得到。

状态空间模型中实际参数估计相关推荐

  1. 状态空间模型中参数的贝叶斯估计

    状态空间模型中参数的贝叶斯估计 (参数估计所有内容) 对于模型中的未知参数 θ∈Rd\bm{\mathbb{\theta \in R^d}}θ∈Rd,贝叶斯方法通常将其建模为先验分布已知的随机变量,且 ...

  2. 卡尔曼滤波 - 状态空间模型中的状态方程

    卡尔曼滤波 - 状态空间模型中的状态方程 flyfish 状态方程和观测方程统称为状态空间模型 位移 位移 = Δ x = x f − x 0 \text { 位移}=\Delta x=x_f-x_0 ...

  3. 【Matlab】状态空间模型的最小化实现 minreal() 函数

    文章目录 含义 例子 Ref 含义 对于单位负反馈,设其开环传函为 GGG,则闭环传函为 Gc=G/(1+G)G_c=G/(1+G)Gc​=G/(1+G),也可以用 feedback() 函数计算,即 ...

  4. matlab中ss函数_matlab状态空间模型(matlab中如何通过ss函数和tf2ss函数将微分方程转化...)...

    matlab中如何通过ss函数和tf2ss函数将微分方程转化... 例如下面的一道题. 方法一: num=[0 0 10 10]; den=[1 6 6 10]; [A,B,C,D]=tf2ss(nu ...

  5. 智能控制中的智能控制系统模型:基于智能模糊状态空间模型的智能控制系统模型

    作者:禅与计算机程序设计艺术 1.简介 在本文中,将阐述智能控制系统模型(Intelligent Control System Model)及其相关理论,包括概率模糊状态空间模型(Probabilis ...

  6. 线性连续时间状态空间模型的离散化及实例

    线性连续时间状态空间模型的离散化(Discretization of Linear Continuous-Time State-Space Models) 1 .状态空间模型 非线性连续时间状态空间模 ...

  7. matlab阶跃响应_状态空间模型及MATLAB指令计算

    一. 基本概念强调 时变控制系统 时变控制系统是指一个或多个系统参数会随着时间变化的系统. 2. 系统状态 系统状态是指表示系统的一组变量,只要知道了这组变量的当前取值情况.知道了输入信号和描述系统动 ...

  8. 多元线性回归模型中多重共线性问题处理方法

    转载自:http://datakung.com/?p=46 多重共线性指自变量问存在线性相关关系,即一个自变量可以用其他一个或几个自变量的线性表达式进行表示.若存在多重共线性,计算自变量的偏回归系数β ...

  9. MATLAB-在命令行估计状态空间模型

    黑箱与结构化状态空间模型估计 黑盒估计 在这种方法中,您可以指定模型顺序,也可以指定配置状态空间矩阵总体结构的其他模型结构属性.您使用数据和模型顺序作为主要输入参数调用ssest.ssregest或n ...

最新文章

  1. 基于 Docker 的 MySQL 导入导出数据
  2. ES6 实用开发技巧
  3. Ardino基础教程 21_LCD1602液晶屏
  4. Java程序员从笨鸟到菜鸟之(七十五)细谈struts2(十四)struts2+ajax实现异步验证...
  5. swig模板 PHP,如何使用nodejs前端模板引擎swig
  6. Go 语言学习总结(5)—— Go 学习笔记总结
  7. php客户端和服务器的值传递
  8. SONiC:为Microsoft全球云提供支持的网络交换机软件
  9. nyoj 366 D的小L(数的全排)
  10. iOS教程:移动终端游戏动画设计的12个原则
  11. 人生的快乐在于自己对生活的态度,快乐是自己的事情
  12. sitf+LK+pnp 识别、跟踪图片,并求三维旋转角度(四) -----LK光流跟踪
  13. Pytorch:手撕ResNet34实现汽车分类
  14. 用STM32CubeMX配置输出PWM信号控制多路舵机(HAL)
  15. 修改Windows Server2003/SQL Server2005的密码后速达软件的配置
  16. 易排通用规划平台,以Excel作为数据源的调用方法与数据文件说明
  17. 最全的android图片加密
  18. C# 制作(PictureBox)圆形头像框 从数据库中读取头像
  19. C51火灾报警器和光源追踪板联调
  20. vue js IOS H5focus无法自动弹出键盘的解决方法

热门文章

  1. python----动态规划
  2. Python paramiko模块基本使用(一)
  3. 字体乱码的时候,可以使用英文下的写法
  4. ADB工具 获取ROOT权限及复制文件方法
  5. RTT的线程同步篇——总结
  6. 计算机基础是高校必修课,高校计算机基础教育教学方式改革
  7. 数据库-MySQL-JDBC-结果集
  8. origin画图_3分钟浏览,Origin绘图中的12个经典问题集锦,早看早知道,躲坑没烦恼!!!...
  9. 7天4场直播,涵盖DBA职业发展必备软实力、Oracle和MySQL技术等
  10. 直播 | 循序渐进 - DM8 数据存储管理