如何使用 R 中的 DSE 包模拟卡尔曼滤波器的后验滤波估计

Posted 2023-02-21

【中文标题】如何使用 R 中的 DSE 包模拟卡尔曼滤波器的后验滤波估计

【英文标题】：How to simulate the posterior filtered estimates of a Kalman Filter using the DSE package in R

【发布时间】：2017-09-18 18:36:32

【问题描述】：

如何使用 DSE 包从 R 中的卡尔曼滤波器模拟中调用后验（精细）状态估计？

我在下面添加了一个示例。假设我创建了一个简单的随机游走状态空间，错误是标准正态分布。该模型是使用 SS 函数创建的，初始化状态和协方差估计为零。因此，理论模型形式为： X(t) = X(t-1) + e(t)~N(0,1) 用于状态演化 Y(t) = X(t) + w(t)~N(0,1)

我们现在按照 Journal of Statistical Software 中“R 中的卡尔曼滤波”文章第 6 页和第 7 页的说明在 R 中实现此功能。首先，我们使用 SS() 函数创建状态空间模型并将其存储在名为 kalman.filter 的变量中：

kalman.filter=dse::SS(F = matrix(1,1,1), 
                  Q = matrix(1,1,1),
                  H = matrix(1,1,1),
                  R = matrix(1,1,1),
                  z0 = matrix(0,1,1),
                  P0 = matrix(0,1,1)
                  )

然后我们使用simulate()从模型表单中模拟100个观察值，并将它们放入一个名为simulate.kalman.filter的变量中：

simulate.kalman.filter=simulate(kalman.filter, start = 1, freq = 1, sampleT = 100)

然后我们使用 l() 对测量结果运行卡尔曼滤波器并将其存储在名为 test 的变量下：

test=l(kalman.filter, simulate.kalman.filter)

从输出中，哪些是我过滤后的估计？

【参考方案1】：

我已经找到了这个问题的答案。

首先，模型的过滤估计值没有在 l() 函数中给出。此功能仅提供领先一步的预测。我的问题的上述框架编码为：

kalman.filter=dse::SS(F = matrix(1,1,1), 
              Q = matrix(1,1,1),
              H = matrix(1,1,1),
              R = matrix(1,1,1),
              z0 = matrix(0,1,1),
              P0 = matrix(0,1,1)
              )
simulate.kalman.filter=simulate(kalman.filter, start = 1, freq = 1, sampleT = 100)
test=l(kalman.filter, simulate.kalman.filter)

领先一步的预测由下式给出：

predictions = test$estimates$pred

一种快速可视化的方法如下：

tfplot(test)

这使您可以根据实际数据快速绘制领先一步的预测。要获得过滤后的估计值，您需要使用同一 dse 包中的 smoother() 函数。它输入状态模型以及数据，在本例中分别为 kalman.filter 和 simulation.kalman.filter。输出是所有时间点的平滑估计。但请注意，它是在考虑完整数据集之后执行此操作的，因此不会在每次观察时都执行此操作。请参见下面的代码。代码的第一行给出了您的平滑估计值，以下几行绘制了示例：

smooth = smoother(test, simulate.kalman.filter)
plot(test$estimates$pred, ylim=c(max(test$estimates$pred,smooth$filter$track,simulate.kalman.filter$outpu), min(test$estimates$pred,smooth$filter$track,simulate.kalman.filter$output)))
points(smooth$smooth$state, col = 3)
points(simulate.kalman.filter$output, col = 4)

上图将您的所有实际数据、模型估计和平滑估计相互绘制。