机器算法验证 - 在 R 中建模多元时间序列计数数据 - 吾爱随笔录

在 R 中建模多元时间序列计数数据

机器算法验证 r 时间序列混合模式计数数据泊松回归

2022-04-13 02:33:21

我正在处理同一类别产品的销售预测（所有属于类别的产品）。数据是具有每日频率的时间序列。经过一番检查，我发现像 Arima 或 Arimax 这样的时间序列模型不太适合我的数据，所以我想尝试一种新方法。 $wine$ $Y_{product_i} \in {0,1,2,..}$

我正在考虑，有计数数据，泊松回归模型，对于单个产品来说，它看起来像这样：

$Y_t|F_{t-1} \sim Poisson(\lambda_{t})$

$\lambda_t = \beta_0 + \beta_1*Y_{t-1} + \alpha_1*\lambda_{t-1} + \eta*X_t$

其中 = t 天的销售额， = t 天的外生回归变量。 $Y_t$ $X_t$

我的第一个问题是：

中有没有实现这个的包？够用吗？ $R$ $glm$

此外，一些产品的时间序列非常稀疏，而另一些产品则至少销售几十年的产品：在第二种情况下考虑对数线性模型，在第一种情况下考虑跨栏模型是否合适？

当时我正在考虑这样一个事实，即我可能需要一个更复杂的模型来模拟不同产品销售之间的协方差结构。不幸的是，在这种情况下，我没有太多的理论和 Rpackages 背景。我正在考虑混合效应模型，但我只是想知道......基本上我想要一个泊松回归模型，我的每个产品都有一个自回归项，但我不知道如何在产品上插入协方差结构。我应该去贝叶斯吗？

$Y_{i_t} \sim Poisson(\lambda_{i_t})$

$\lambda_{i_t} = \beta_{0_i} + \beta_{i_1}*Y_{i_{t-1}} +\alpha_{1_i}*\lambda_{i_{t-1}} + \eta_i*X_{i_t}$

$i=1,..,numberofproducts$

所以我的第二个问题是关于找到关于多元时间序列计数数据以及 R 中的包的理论观点。

1个回答

我找到了回答我的问题的参考：https ://arxiv.org/pdf/1405.3738.pdf 。

模型相当复杂，这里是状态空间表示：

所以，假设我在 1,..,T 个时间段内研究了 L 种不同的产品。

$Y_{l,t} \sim z*\delta_0 + (1-z)NB(exp(\widetilde{\eta}_{l,t}),alpha_l)$ 是产品 l 在时间 t 的分布

$\widetilde{\eta}_{l,t} = \eta_{l,t} + X_{l,t}\theta_l$ 这是产品 l 在时间 t 的销售平均值的对数，保证它是正的。

$\eta_{l,t} = \mu_l + \phi_l(\eta_{l,t-1}-\mu) + \epsilon_{l,t}$

$\epsilon_{l,t} \sim N(0,\frac{1}{\tau_l})$

其他先验和超先验在下图中：

PS 现在我正在尝试编写 JAGS 代码，任何帮助将不胜感激！（https://stackoverflow.com/questions/40528715/runtime-error-in-jags）

编辑：

这是 JAGS 代码：

model{


#hyperpriors 4
alpha_star ~ dunif(0.001,0.1)
tau_mu_star ~ dunif(1,10)
mu_star ~ dnorm(0,0.5)
beta_tau ~ dunif(2,25)
beta_0_tau ~ dunif(1,10)
beta_theta ~ dunif(2,25)
phiminus ~ dunif(1,50)
k_tau ~ dunif(5,10)
k_0_tau ~ dunif(1,5)
pointmass_0 ~ dnorm(0,10000)
k_theta ~ dunif(5,10)
phiplus ~ dunif(1,600)
theta_star ~ dmnorm(b0,B0)
#17
for (l in 1:L){
z[l] ~ dbeta(0.5,0.5)
phi[l] ~ dbeta(phiplus + phiminus, phiminus)
tau[l] ~ dgamma(k_tau,beta_tau)
tau_theta[l] ~ dgamma(k_tau,beta_tau)
mu[l] ~ dnorm(mu_star, tau_mu_star)
alpha[l] ~ dexp(alpha_star)
eps[1,l] ~ dnorm(0,tau[l])
eta[1,l] = mu_star + eps[1,l]
theta[l,1:8] ~ dmnorm(theta_star,thetavar*tau_theta[l])
#y[1,l] ~ inprod(1-z[l],dnegbin(exp(eta[1,l]),alpha[l]))
y[1,l] ~ dnegbin(exp(eta[1,l]),alpha[l])

#y[1,l] ~ dnegbin(exp(eta[1,l]),alpha[l])
ystar[1,l] ~ dnorm(z[l]*pointmass_0 + inprod((1-z[l]),y[1,l]),100000)
}

for (i in 2:N){

for (l in 1:L){
eps[i,l] ~ dnorm(0,tau[l])
}

    for(l in 1:L){
    eta[i,l] = mu[l]+ phi[l]*(eta[i-1,l]-mu[l]) + eps[i,l] 
    eta_star[i,l]= eta[i,l] + inprod(c(x[i,l],xshared[i,]),t(theta[l,]))
    #observations
#kobe[i,l] ~ dnegbin(dexp(eta_star[i,l]),alpha[l])
#   #y[i,l] = inprod(1-z[l],kobe[i,l])
    #y[i,l] ~ inprod(1-z[l],dnegbin(exp(eta_star[i,l]),alpha[l]))
    #y[i,l] ~ dnegbin(exp(eta_star[i,l]),alpha[l])
    y[i,l] ~ dnegbin(exp(eta_star[i,l]),alpha[l])
    ystar[i,l] ~ dnorm(z[l]*pointmass_0 + inprod((1-z[l]),y[i,l]),100000)
}
}

}

我使用 runjags 从 R 调用：

 parsamples <- run.jags('jags_model.txt', data=forJags, monitor=c('y','theta'), sample=100, method='rjparallel')

其它你可能感兴趣的问题

上一篇马尔可夫随机场与隐马尔可夫模型下一篇定义先验多项式回归。`MCMCglmm` 的案例研究