使用解释变量对时间序列数据进行建模
Posted
技术标签:
【中文标题】使用解释变量对时间序列数据进行建模【英文标题】：Modelling time series data with explanatory variable 【发布时间】：2019-11-03 15:22:31 【问题描述】：我有以下数据，其中
y 是一些产品的篮子的平均价格 x 是汇率。 period y x
1 201501 1530 2.49
2 201502 1450 2.62
3 201503 1637 2.77
4 201504 1404 2.84
5 201505 1442 2.82
6 201506 1442 2.89
7 201507 1518 2.88
8 201508 1492 3.05
9 201509 1743 3.21
10 201510 1902 3.14
11 201511 1855 3.07
12 201512 1879 3.12
13 201601 2018 3.21
14 201602 2117 3.15
15 201603 2002 3.09
16 201604 1837 3.04
17 201605 1902 3.14
18 201606 1910 3.12
19 201607 2162 3.16
20 201608 2183 3.17
21 201609 2100 3.17
22 201610 2122 3.28
23 201611 2461 3.51
24 201612 2250 3.73
25 201701 2466 4.00
26 201702 2212 3.93
27 201703 2424 3.93
28 201704 2477 3.91
29 201705 2402 3.82
30 201706 2360 3.77
31 201707 2475 3.81
32 201708 2690 3.76
33 201709 2655 3.70
34 201710 2889 3.92
35 201711 2683 4.15
36 201712 2674 4.12
37 201801 2695 4.03
38 201802 2707 4.04
39 201803 2728 4.15
40 201804 2607 4.33
41 201805 2917 4.71
42 201806 2946 4.94
43 201807 3031 5.08
44 201808 3224 6.20
45 201809 3962 6.76
46 201810 4043 6.25
47 201811 3805 5.76
48 201812 3607 5.67
49 201901 3694 5.74
50 201902 3566 5.63
51 201903 3541 5.83
52 201904 3350 6.15
所以为了揭示 y 和 x（和 x 的滞后）之间的关系，我做了以下分析：
library(DataCombine)
#data
data2<-structure(c(201501, 201502, 201503, 201504, 201505, 201506, 201507,
201508, 201509, 201510, 201511, 201512, 201601, 201602, 201603,
201604, 201605, 201606, 201607, 201608, 201609, 201610, 201611,
201612, 201701, 201702, 201703, 201704, 201705, 201706, 201707,
201708, 201709, 201710, 201711, 201712, 201801, 201802, 201803,
201804, 201805, 201806, 201807, 201808, 201809, 201810, 201811,
201812, 201901, 201902, 201903, 201904, 1530, 1450, 1637, 1404,
1442, 1442, 1518, 1492, 1743, 1902, 1855, 1879, 2018, 2117, 2002,
1837, 1902, 1910, 2162, 2183, 2100, 2122, 2461, 2250, 2466, 2212,
2424, 2477, 2402, 2360, 2475, 2690, 2655, 2889, 2683, 2674, 2695,
2707, 2728, 2607, 2917, 2946, 3031, 3224, 3962, 4043, 3805, 3607,
3694, 3566, 3541, 3350, 2.49, 2.62, 2.77, 2.84, 2.82, 2.89, 2.88,
3.05, 3.21, 3.14, 3.07, 3.12, 3.21, 3.15, 3.09, 3.04, 3.14, 3.12,
3.16, 3.17, 3.17, 3.28, 3.51, 3.73, 4, 3.93, 3.93, 3.91, 3.82,
3.77, 3.81, 3.76, 3.7, 3.92, 4.15, 4.12, 4.03, 4.04, 4.15, 4.33,
4.71, 4.94, 5.08, 6.2, 6.76, 6.25, 5.76, 5.67, 5.74, 5.63, 5.83,
6.15), .Dim = c(52L, 3L), .Dimnames = list(NULL, c("period",
"y", "x")), .Tsp = c(1, 5.25, 12),
class = c("mts", "ts", "matrix"))
# deaseasonal data using loess procedure
# model assumed to be multiplicative ->
# so seasonal coefficients obtained after taking logarithm
data2<-ts(data1, frequency = 12)
lprod<-log(data2[,2])
decomp<-stl(lprod, s.window="periodic")
decomp<-decomp$time.series
season<-exp(decomp[,1])
trend<-exp(decomp[,2])
rand<-exp(decomp[,3])
#deasonal y value
desdata<-trend*rand
#obtaining lags(1 and 2) of explonatary variable x
ex_var<-as.data.frame(data2)
ex_var<-slide(ex_var, Var='x', NewVar = "x1", slideBy = -1)
ex_var<-slide(ex_var, "x", NewVar = "x2", slideBy = -2)
ex_var<-slide(ex_var, "x", NewVar = "x3", slideBy = -2)
#delete firts two rows
ex_var<-ex_var[-c(1:2),]
#regression
#I also include x variable at time t. My aim is just to obtaion the relation
#After some trials, the below models is fitted
myreg<-lm(formula=y~-1+x+x1, data=ex_var)
#fitted values and deviations
fitted<-myreg$fitted.values
dev<-(fitted/ex_var$y-1)*100
所以，我只是将 y 去季节化并在 x 和 x1 上回归 y。
然后，我得到了拟合值和偏差。
dev
3 1.8692282
4 24.5705019
5 23.3397058
6 23.4631450
7 19.3359159
8 22.9155774
9 11.2428038
10 5.2958604
11 5.5962845
12 2.9269927
13 -2.2926331
14 -5.3226638
15 -1.7643976
16 5.0966933
17 1.1112870
18 2.9722925
19 -9.1685959
20 -9.1115852
21 -5.2963860
22 -5.4530041
23 -14.8963143
24 -0.5730753
25 -3.3622134
26 12.9454251
27 1.7153551
28 -0.5895710
29 1.5281120
30 1.2122606
31 -4.1790659
32 -11.4373380
33 -11.5113009
34 -18.4386679
35 -6.9727496
36 -2.8112683
37 -4.6213799
38 -6.5419446
39 -6.4478475
40 0.9688325
41 -4.7986412
42 1.5460511
43 2.9863159
44 4.3845503
45 0.4257484
46 2.8904249
47 1.0008665
48 -0.2121615
49 -3.4013337
50 0.4939980
51 0.6482636
52 10.7024816
似乎偏差是高度自相关的。
所以，我想我必须构建一个包含 y 滞后的模型。但我不确定。我需要取对数的差异吗？还是我应该遵循其他方式？
我不是时间序列数据方面的专家。我真的被这些数据困住了。我会很高兴得到任何帮助。非常感谢您从现在开始节省您的时间。
【问题讨论】：
【参考方案1】：我无法检测到季节性。
对于有滞后的线性回归，您可以使用dynlm，但您还应该考虑使用 (pseudo) ARIMAX 模型。
library(zoo)
library(dynlm)
# Convert to a multivariate zoo object
ym <- as.yearmon(as.character(data2[,"period"]), "%Y%m")
tt.zoo <- zoo(data2[,c("y", "x")], order.by=ym)
# No significant periodicity
tt.d <- diff(tt.zoo)
ccf(tt.d[,"x"], tt.d[,"y"])
acf(tt.d[,"y"])
acf(tt.d[,"x"])
# Dynamic Linear Models of various complexity
dlm0 <- dynlm(y ~ L(y, 1),
data=tt.zoo, start=start(tt.zoo)+1/12, end=end(tt.zoo))
dlm1 <- dynlm(y ~ x,
data=tt.zoo, start=start(tt.zoo)+1/12, end=end(tt.zoo))
dlm2 <- dynlm(y ~ x + L(x, 1),
data=tt.zoo, start=start(tt.zoo)+1/12, end=end(tt.zoo))
dlm3 <- dynlm(y ~ x + L(x, 1) + L(y, 1),
data=tt.zoo, start=start(tt.zoo)+1/12, end=end(tt.zoo))
# Residuals look reasonably OK. Maybe a slight curve?
plot(residuals(dlm0), col=2)
lines(residuals(dlm1), col=3)
lines(residuals(dlm2), col=4)
lines(residuals(dlm3), col=5)
# Fits are pretty OK
plot(tt.zoo[,"y"], col="grey", lwd=2)
lines(fitted(dlm0), col=2)
lines(fitted(dlm1), col=3)
lines(fitted(dlm2), col=4)
lines(fitted(dlm3), col=5)
# dlm1 and dlm2 have significant acf, due to y's non-stationarity
acf(residuals(dlm0))
acf(residuals(dlm1))
acf(residuals(dlm2))
acf(residuals(dlm3))
# dlm3 seems like it makes good use of the extra marameters
anova(dlm0, dlm1, dlm2, dlm3)
AIC(dlm0, dlm1, dlm2, dlm3)
# df AIC
# dlm0 3 678.6580
# dlm1 3 699.1142
# dlm2 4 686.8249
# dlm3 5 660.2720
# Altermative 'ARIMAX' model with lagged external regressor
tt.l <- cbind(tt.zoo, xl=lag(tt.zoo[,2], -1))
tt.l <- tt.l[complete.cases(tt.l),]
ax1 <- arima(tt.l[,1], order=c(0, 1, 0), xreg=tt.l[,2])
ax2 <- arima(tt.l[,1], order=c(0, 1, 0), xreg=tt.l[,2:3])
ax3 <- arima(tt.l[,1], order=c(1, 1, 0), xreg=tt.l[,2:3])
# Again the extra parameters seems to be justified
AIC(ax1, ax2, ax3)
# df AIC
# ax1 2 656.1224
# ax2 3 649.2153
# ax3 4 645.8715
# Checking that the residuals are sufficiently white noise-like
library(forecast)
checkresiduals(ax1)
checkresiduals(ax2)
checkresiduals(ax3)
【讨论】：
非常感谢您详细的退货。我猜，在 arimax 部分，tt.l 处的 x1 不是我猜的 x 的滞后值。它似乎是 x 的前向值。非常感谢。 @oercim:tt.l 有三个变量：y 和 x 与您的数据相同，xl 是 x，但落后一位。比较head(tt.zoo); head(tt.l)，你会看到它是如何排列的。这意味着ax2 和ax3 使用未滞后和滞后x 作为外部回归器，而ax1 仅使用x。
@AksalA，当我查看时，tt.l 的第一行 - 2015 年 1 月，xl 的值为 2.62，它是 2015 年 2 月的值。它似乎是一个前进（或领先的价值）。我猜第一行的 xl 应该是空的，因为我们不知道 2014 年 12 月的值。
啊，是的，你是对的。 lag() 落后于索引，而不是数据。我会改正的。
，非常感谢。这对我来说是一个非常有用和有帮助的答案。最好的问候。以上是关于使用解释变量对时间序列数据进行建模的主要内容，如果未能解决你的问题，请参考以下文章
Python对商店数据进行lstm和xgboost销售量时间序列建模预测分析|附代码数据
Python对商店数据进行lstm和xgboost销售量时间序列建模预测分析|附代码数据