Python ARIMA 外生变量样本外分析

11 投票
3 回答
26035 浏览
提问于 2025-04-18 15:23

我正在尝试使用Python的statsmodels库中的ARIMA包来预测一个时间序列,并且想要加入一个外部变量,但我不知道在预测步骤中该如何正确地插入这个外部变量。你可以在这里查看相关文档。

import numpy as np
from scipy import stats
import pandas as pd

import statsmodels.api as sm

vals = np.random.rand(13)
ts = pd.TimeSeries(vals)
df = pd.DataFrame(ts, columns=["test"])
df.index = pd.Index(pd.date_range("2011/01/01", periods = len(vals), freq = 'Q'))

fit1 = sm.tsa.ARIMA(df, (1,0,0)).fit()
#this works fine:
pred1 = fit1.predict(start=12, end = 16)
print(pred1)

Out[32]: 
2014-03-31    0.589121
2014-06-30    0.747575
2014-09-30    0.631322
2014-12-31    0.654858
2015-03-31    0.650093
Freq: Q-DEC, dtype: float64

现在加入一个趋势的外部变量。

exogx = np.array(range(1,14))
#to make this easy, let's look at the ols of the trend (arima(0,0,0))
fit2 = sm.tsa.ARIMA(df, (0,0,0),exog = exogx).fit()
print(fit2.params)

const    0.555226
x1       0.013132
dtype: float64

print(fit2.fittedvalues)

2011-03-31    0.568358
2011-06-30    0.581490
2011-09-30    0.594622
2011-12-31    0.607754
2012-03-31    0.620886
2012-06-30    0.634018
2012-09-30    0.647150
2012-12-31    0.660282
2013-03-31    0.673414
2013-06-30    0.686546
2013-09-30    0.699678
2013-12-31    0.712810
2014-03-31    0.725942
Freq: Q-DEC, dtype: float64

注意,正如我们所预期的,这是一条趋势线,随着时间的推移,每增加一个单位就增加0.013132(当然,这些数据是随机的,所以如果你运行它,数值会不同,但趋势的方向是一样的)。因此,下一个值(当时间=14时)应该是0.555226 + 0.013132*14 = 0.739074。

#out of sample exog should be (14,15,16)
pred2 = fit2.predict(start = 12, end = 16, exog = np.array(range(13,17)))
print(pred2)
2014-03-31    0.725942
2014-06-30    0.568358
2014-09-30    0.581490
2014-12-31    0.594622
2015-03-31    0.765338
Freq: Q-DEC, dtype: float64

所以,2014年3月31日的预测是正确的(这是最后一个样本),但2014年6月30日又回到了开始(t = 1),不过注意2015年3月31日(实际上,总是预测的最后一个观察值,无论预测的时间范围)是t = 16(也就是说,(值 - 截距)/beta = (0.765338 - 0.555226)/0.013132)。

为了让这个更清楚,注意当我增加x矩阵的值时会发生什么。

fit2.predict(start = 12, end = 16, exog = np.array(range(13,17))*10000)
Out[41]: 
2014-03-31       0.725942
2014-06-30       0.568358
2014-09-30       0.581490
2014-12-31       0.594622
2015-03-31    2101.680532
Freq: Q-DEC, dtype: float64

你看到2015年3月31日的值飙升了,但其他的xmat值却没有被考虑到?我到底哪里做错了???

我尝试了我知道的所有方法来传递外部变量(改变维度、把外部变量变成矩阵、让外部变量的长度等于输入加上预测范围等等)。如果有任何建议,我将非常感激。

我使用的是Anaconda2.1的2.7版本,numpy 1.8.1,scipy 0.14.0,pandas 0.14.0,statsmodels 0.5.0。

我在64位的Windows 7和64位的CentOS上验证了这个问题。

另外,有几点需要说明。我使用ARIMA是为了ARIMA的功能,上面的内容只是为了说明(也就是说,我不能“仅仅使用OLS...”,我想这会被建议)。由于项目的限制,我也不能“仅仅使用R”(更一般来说,Spark的基础版本对R的支持不足)。

这里是代码中有趣的部分,如果你想自己尝试一下。

import numpy as np
from scipy import stats
import pandas as pd
import statsmodels.api as sm

vals = np.random.rand(13)
ts = pd.TimeSeries(vals)
df = pd.DataFrame(ts, columns=["test"])
df.index = pd.Index(pd.date_range("2011/01/01", periods = len(vals), freq = 'Q'))

exogx = np.array(range(1,14))
fit2 = sm.tsa.ARIMA(df, (0,0,0),exog = exogx).fit()
print(fit2.fittedvalues)
pred2 = fit2.predict(start = 12, end = 16, exog = np.array(range(13,17))*10000)
print(pred2)

3 个回答

1

如果有人在使用 .forecast() 这个方法,我用它做了一步预测,效果不错。

  • history 是用来训练的数组

  • exog 是外部变量的数组

  • Y_exog_test 是对应的外部变量,用于测试的样本。

把它换成 .SARIMAX(),应该就能正常工作了。

model = sm.tsa.statespace.ARIMA(history, trend='c', order=(1,1,1),seasonal_order=(0,1,0,24),exog=yexog)
model_fit = model.fit()
predicted = model_fit.forecast(step=1,exog=[[Y_exog_test]], dynamic=True)
2

在调整 fit2 的时候,你已经提到过外部变量,所以不需要再重复了:

exogx = np.array(range(1,5)) # I think you will need 4 exegeneous variables to perform an ARIMAX(0,0,0) since you want out of sample forecast with 4 steps ahead
fit2 = sm.tsa.ARIMA(df, (0,0,0),exog = exogx).fit()
# if you want to do an out-of-sample-forecast use fit2.forecast(steps) instead
#I would do this
pred = fit2.forecast(steps = 4)
fcst_index = pd.date_range(start = df.shift(1,'10T').index[-1]  , periods = 4, freq = '10T')
fcst_serie = pd.Series(data = pred1[0], index = fcst_index)
print fcst_serie

希望这能帮到你!这篇文章很棒。我之前从来没有在 ARIMA 中尝试过外部变量,但有些论文说无论你在哪个领域使用,它其实都不是很相关(如果需要的话我可以帮你找这些论文,或者你也可以自己去谷歌一下)。

4

这个问题最好发到github的问题追踪器上。我已经提交了一个工单

如果能在那儿提交工单会更好,不然我可能会忘记这件事。这段时间我挺忙的。

在k_ar等于0的特殊情况下,逻辑上有个bug。现在应该修好了。如果你能试试这个修复,告诉我结果如何。如果不能的话,我可以做一些更严格的测试,然后合并这个修复。

在spark上使用Statsmodels?我很感兴趣。

撰写回答