R语言状态空间模型和卡尔曼滤波预测酒精死亡人数时间序列

您所在的位置:网站首页 状态空间观测器 R语言状态空间模型和卡尔曼滤波预测酒精死亡人数时间序列

R语言状态空间模型和卡尔曼滤波预测酒精死亡人数时间序列

2023-03-24 03:01| 来源: 网络整理| 查看: 265

原文链接:http://tecdat.cn/?p=22665 摘要

状态空间建模是一种高效、灵活的方法,用于对大量的时间序列和其他数据进行统计推断。本文介绍了状态空间建模,其观测值来自指数族,即高斯、泊松、二项、负二项和伽马分布。在介绍了高斯和非高斯状态空间模型的基本理论后,提供了一个泊松时间序列预测的说明性例子。最后,介绍了与拟合非高斯时间序列建模的其他方法的比较。

绪论

状态空间模型为几种类型的时间序列和其他数据的建模提供了一个统一的框架。结构性时间序列、自回归综合移动平均模型(ARIMA)、简单回归、广义线性混合模型和三次样条平滑模型只是一些可以表示为状态空间模型的统计模型的例子。最简单的一类状态空间模型是线性高斯状态空间模型(也被称为动态线性模型),经常被用于许多科学领域。

高斯状态空间模型

本节将介绍有关高斯状态空间模型理论的关键概念。由于卡尔曼滤波(Kalman filtering)背后的算法主要是基于Durbin和Koopman(2012)以及同一作者的相关文章。对于具有连续状态和离散时间间隔的线性高斯状态空间模型t=1, . . . ,n,我们有

7469b6bc9c6390fb4e91b317fbd3c70e.png

其中t∼N(0,Ht),ηt∼N(0,Qt)和α1∼N(a1,P1)相互独立。我们假设yt是一个p×1,αt+1是一个m×1,ηt是一个k×1的向量。α = (α > 1 , . . . , α> n ) >,同样y = (y > 1 , . . , y> n ) >。

状态空间建模的主要目标是在给定观测值y的情况下获得潜状态α的知识。这可以通过两个递归算法实现,即卡尔曼滤波和平滑算法。从卡尔曼滤波算法中,我们可以得到先行一步的预测结果和预测误差

5afd8b8c93cbe9ce33c5072f6b11d56a.png

和相关的协方差矩阵

16a98e085322fa93e95fc886705a863c.png

利用卡尔曼滤波的结果,我们建立了状态平滑方程,在时间上向后运行,产生了

d40b8cbf9fb42d61e173c05cac859583.png

对于干扰项t和ηt,对于信号θt = Ztαt,也可以计算类似的平滑估计。

高斯状态空间模型的例子

现在通过例子来说明。我们的时间序列包括1969-2007年40-49岁年龄组每年每10万人中酒精相关的死亡人数(图1)。数据取自统计局。对于观测值 y1, ... . , yn,我们假设在所有t = 1, . . . , n,其中μt是一个随机游走的漂移过程

2cfcd31d8b74987e6a03ecea4677f30a.png

ηt∼N(0, σ2 η)。假设我们没有关于初始状态μ1或斜率ν的先验信息。这个模型可以用状态空间的形式来写,定义为

f4f91ce8d0bf410cd4f0a60ec22b83b7.png

在KFAS中,这个模型可以用以下代码来写。为了说明问题,我们手动定义所有的系统矩阵,而不采用默认值。

R> Zt  model_gaussian  Model(deaths ~ -1 + + distribution = "poisson")

与之前的高斯模型相比,我们现在需要用参数distribution(默认为 "高斯")来定义观测数据的分布。我们还通过参数u来定义暴露项,并使用a1和P1的默认值。在这个模型中,只有一个未知参数,即σ 2 η。这个参数被估计为0.0053,但是高斯模型和泊松模型之间σ 2 η的实际值不能直接比较,因为不同模型对µt的解释不同。泊松模型的斜率项估计为0.022,标准误差为1.4×10-4,对应于死亡人数每年增加2.3%。

图2显示了以高斯过程(蓝色)和泊松过程(红色)为模型(每10万人的死亡人数)的平滑估计。

df0b6dd50b1ae74a318b47db906b8527.png

任意的状态空间模型

通过结合前面的方法,可以相对容易地构建大量的模型。对于这样做还不够的情况,可以通过直接定义系统矩阵来构建任意状态空间模型。作为一个例子,我们修改了之前的泊松模型,增加了一个额外的白噪声项,试图捕捉数据的可能的过度离散。现在我们的泊松强度模型是ut exp(µt + t),即

23ae0c328aeb70be0cbc6feef14662fc.png

其中ηt∼N(0, σ2 η)如前,t∼N(0, σ2)。这个模型可以用状态空间的形式来写,定义为

Model(deaths ~ trend(2, Q = list(NA, 0)) + distribution = "poisson")

由于模型包含P1中的未知参数,我们需要提供一个特定的模型更新函数。

R> update  updatefn  varcordel\["Q",  "custom"\]

25255e31c5fd7d7b2d8eb0c017a7726c.png

b826a834619dda2d7a90bfb1a39c390b.png

状态空间模型的参数估计通常工作量很大,因为似然面包含多个最大值,从而使优化问题高度依赖于初始值。通常情况下,未知参数与未观察到的潜在状态有关,如本例中的协方差矩阵,几乎没有先验知识。

因此,要猜出好的初始值是很有挑战性的,特别是在更复杂的环境中。因此,在可以合理地确定找到适当的最优值之前,建议使用多种初始值配置,可能有几种不同类型的优化方法。这里我们使用观察到的系列的协方差矩阵作为协方差结构的初始值。在非高斯模型的情况下,另一个问题是,似然计算是基于迭代程序的,该程序使用一些终止条件(如对数似然的相对变化)停止,因此对数似然函数实际上包含一些噪声。这反过来又会影响BFGS等方法的梯度计算,在理论上可以得到不可靠的结果。因此,有时建议使用无导数的方法,如Nelder-Mead。另一方面,BFGS通常比Nelder-Mead快得多,因此我更愿意先尝试BFGS,至少在初步分析中。我们可以计算出状态的平滑估计。

R> out 


【本文地址】


今日新闻


推荐新闻


CopyRight 2018-2019 办公设备维修网 版权所有 豫ICP备15022753号-3