0 SARIMAX模型时间序列分析步骤

1. 用pandas处理时序数据

2. 检验时序数据的平稳性

3. 将时序数据平稳化

4. 确定order 的 p.d.q值

5. 确定season_order的四个值

6. 应用SARIMAX模型对时序数据进行预测

其实SARIMAX比ARIMA模型就多了个season_order参数的确定，但也是这里最费时间的一个步骤

1 将数据转化成为时序数据

先一股脑导入一下工具包

import pandas as pd
import datetime
import matplotlib.pyplot as plt
from pylab import mpl
mpl.rcParams['font.sans-serif']=['SimHei']
import seaborn as sns
import statsmodels.tsa.stattools as ts
import statsmodels.api as sm
from statsmodels.tsa.arima_model import ARIMA
from statsmodels.stats.diagnostic import unitroot_adf
from statsmodels.graphics.tsaplots import plot_acf, plot_pacf
import itertools
import warnings
import numpy as np
from statsmodels.tsa.seasonal import seasonal_decompose

#读取数据
data = pd.read_csv('factor.csv')
data.index = pd.to_datetime(data['date'])
data.drop(['date'], axis=1, inplace=True)
data = data.result
data.head()

#数据大致情况展示
data.plot(figsize=(12,8))
plt.legend(bbox_to_anchor=(1.25, 0.5))
plt.title('result')
sns.despine()
plt.show()

2 序列平稳性检测

#数据平稳性检测 因为只有平稳数据才能做时间序列分析
def judge_stationarity(data_sanya_one):dftest = ts.adfuller(data_sanya_one)print(dftest)dfoutput = pd.Series(dftest[0:4], index=['Test Statistic','p-value','#Lags Used','Number of Observations Used'])stationarity = 1for key, value in dftest[4].items():dfoutput['Critical Value (%s)'%key] = value if dftest[0] > value:stationarity = 0print(dfoutput)print("是否平稳(1/0): %d" %(stationarity))return stationarity
stationarity = judge_stationarity(data)

3 序列平稳化

#若不平稳进行一阶差分
if stationarity == 0:data_diff = data.diff()data_diff = data_diff.dropna()plt.figure()plt.plot(data_diff)plt.title('一阶差分')plt.show()#再次进行平稳性检测
stationarity = judge_stationarity(data_diff)

4 做一下季节性分解看看有没有季节性

#季节性分解
decomposition = seasonal_decompose(data,freq=28)
trend = decomposition.trend
seasonal = decomposition.seasonal
residual = decomposition.residplt.figure(figsize=[15, 7])
decomposition.plot()
print("test: p={}".format(ts.adfuller(seasonal)[1]))#季节平稳性检测
stationarity = judge_stationarity(residual)

5 对order参数p、q定阶，下面两种都可以，但图我还是看不来，下一种傻瓜式（稍微费时间）

#画ACF图和PACF图来确定p、q值
from statsmodels.graphics.tsaplots import plot_acf, plot_pacfdef draw_acf_pacf(ts,lags):f = plt.figure(facecolor='white')ax1 = f.add_subplot(211)plot_acf(ts,ax=ax1,lags=lags)  #lags 表示滞后的阶数，值为30，显示30阶的图像ax2 = f.add_subplot(212)plot_pacf(ts,ax=ax2,lags=lags)  plt.subplots_adjust(hspace=0.5)plt.show()
draw_acf_pacf(myts_diff,30)

#对模型p,q进行定阶
warnings.filterwarnings("ignore") # specify to ignore warning messages
from statsmodels.tsa.arima_model import ARIMA pmax = int(5)    #一般阶数不超过 length /10
qmax = int(5)
bic_matrix = []
for p in range(pmax +1):temp= []for q in range(qmax+1):try:temp.append(ARIMA(data, (p, 1, q)).fit().bic)except:temp.append(None)bic_matrix.append(temp)bic_matrix = pd.DataFrame(bic_matrix)   #将其转换成Dataframe 数据结构
p,q = bic_matrix.stack().idxmin()   #先使用stack 展平， 然后使用 idxmin 找出最小值的位置
print(u'BIC 最小的p值 和 q 值：%s,%s' %(p,q))  #  BIC 最小的p值 和 q 值：0,1

6 通过网格搜索对seasonal_order进行定阶

#通过网格搜索对seasonal_order进行定阶,目前就是pdq=011，seasonal_order=2, 2, 1, 52效果比较好，RMSE=202.4582
def get_ARIMA_params(data, pdq, m=12):p = d = q = range(0, 3)seasonal_pdq = [(x[0], x[1], x[2], m) for x in list(itertools.product(p, d, q))]score_aic = 1000000.0warnings.filterwarnings("ignore") # specify to ignore warning messagesfor param_seasonal in seasonal_pdq:mod = sm.tsa.statespace.SARIMAX(data,order=pdq,seasonal_order=param_seasonal,enforce_stationarity=False,enforce_invertibility=False)results = mod.fit()print('x{}12 - AIC:{}'.format(param_seasonal, results.aic))if results.aic < score_aic:score_aic = results.aicparams = param_seasonal, results.aicparam_seasonal, results.aic = paramsprint('x{}12 - AIC:{}'.format(param_seasonal, results.aic))
pdq = [0, 1, 1]
get_ARIMA_params(data, pdq, m=52)

这上面最关键的是这个m值怎么设定，我设置默认为12。m是季节周期，参考别人代码的时候，月度数据的m为12，那我这里以周为单位，季节周期应该是52吧。一年12个月，52个星期，这个逻辑应该没有问题。

7 根据定阶参数进行模型拟合

mod = sm.tsa.statespace.SARIMAX(data,order=(0, 1, 1),seasonal_order=(2, 1, 2, 52),enforce_stationarity=False,enforce_invertibility=False)
results = mod.fit()
print(results.summary().tables[1])
results.plot_diagnostics(figsize=(15, 12))
plt.show()

这里模型拟合的时候用的数据的原始数据，而不是差分后的数据，因为order参数中已经设置了d为1，在拟合的时候会自动进行一阶差分，并在预测的时候对预测结果进行差分还原。

8 对预测值和真实值作图，并计算RMSE值作为评估参数

predict_ts = results.predict(tpy='levels')  #tpy='levels'直接预测值，没有的话预测的是差值
myts = data[predict_ts.index]  # 过滤没有预测的记录predict_ts.plot(color='blue', label='Predict',figsize=(12,8))myts.plot(color='red', label='Original',figsize=(12,8))plt.legend(loc='best')
plt.title('RMSE: %.4f'% np.sqrt(sum((predict_ts-myts)**2)/myts.size))
plt.show()

9 向后对forecast值作图

steps = 20
start_time = myts.index[-1]
forecast_ts = results.forecast(steps)fore = pd.DataFrame()
fore['date'] = pd.date_range(start=start_time ,periods=steps, freq='7D')
fore['result'] = pd.DataFrame(forecast_ts.values)
fore.index = pd.to_datetime(fore['date'])predict_ts['2019/1/18':].plot(color='blue', label='Predict',figsize=(12,8))
myts['2019/1/18':].plot(color='red', label='Original',figsize=(12,8))
fore.result.plot(color='black', label='forecast',figsize=(12,8))plt.legend(loc='best')
plt.show()

这里有个函数pd.date_range是专门用于产生时间序列索引的，start = 开始时间，end = 结束时间，periods=时间索引的个数，freq=‘7D’表示7天为一个时间索引间隔，也可以是'7W'七周，'M'一个月等等。

由于预测的数据没有时间索引，只有序号所以我要在这给他生成时间索引，并合并到dataframe，这样就可以和其他值一起在图像上展示了。

最后forecast的效果还是可以的嘛，保存forecast文件

fore.to_csv('forecast_20steps.csv')

季节性时间序列分析-SARIMAX模型的python实现相关推荐

python时间序列分析航空旅人_时间序列分析-ARIMA模型（python）
时间序列概念:在生产和科学研究中,对某一个或者一组变量进行观察测量,将在一系列时刻所得到的离散数字组成的序列集合,称之为时间序列.时间序列分析是根据系统观察得到的时间序列数据,通过曲线拟合和参数估 ...
时间序列分析ARMA模型原理及Python statsmodels实践（上）
目录 1. 时间序列及相关基本概念 1.1. 时间序列分解 1.2. 时间平稳序列 1.3. 自相关与自相关函数(ACF) 1.4. 白噪声及Ljung-Box检验 1.4.1. 白噪声 1.4.2. ...
arima 数据预处理_时间序列分析|ARIMA模型分步骤解析及R中实践
你是否想要做时间序列分析,但却不知道代码怎么写? 你是否不清楚时间序列分析各种模型该在什么情况下使用? 本文将针对以上两个问题,带你入门时间序列分析~ 等等! 不止'入门' 读完这篇,你立即就能在R中 ...
时间序列分析-ARIMA模型
提示:文章写完后,目录可以自动生成,如何生成可参考右边的帮助文档时间序列分析-ARIMA模型概述一.时间序列平稳性 1.严平稳 2.弱平稳二.建模步骤 1.平稳性检验 2.人工p 和 q 阶数 ...
covariance matrix r语言_时间序列分析|ARIMAX模型分步骤详解和R中实践
这是关于时间序列的第N篇文章,本文将介绍ARIMAX模型,简单来说就是在ARIMA的基础上增加一个外生变量.ARIMAX和ARIMA相比在理论上没有太多新的内容,所以本文直接介绍在R里怎么一步一步跑A ...
python arima模型_时间序列分析 ARIMA模型 Python（2）
最近做了一个时间序列分析的项目.时间序列分析不同于以前的项目--看一下相关的库怎么用,就可以快速上线应用.它是非常需要你的基础知识的,Hamilton关于<时间序列分析>方面的知识,写了厚 ...
时间序列分析之ARIMA上手-Python
概念时间序列时间序列(或称动态数列)是指将同一统计指标的数值按其发生的时间先后顺序排列而成的数列.时间序列分析的主要目的是根据已有的历史数据对未来进行预测. 时间序列分析时间序列分析是根据系统观 ...
数学建模--时间序列分析、模型预测
学习自:b站清风数学建模第1部分_时间序列分析的概念与时间序列分解模型时间序列时间序列的基本概念区分时期和时点序列时期序列适用于灰色预测模型时间序列分解时间要小于1年,才能进行时间序列 ...
SVR，时间序列分析的评价指标，python数据标准化
知识点 SVR 参考支持向量机(SVM)是一种分类算法,但是也可以做回归,根据输入的数据不同可做不同的模型(若输入标签为连续值则做回归,若输入标签为分类值则用SVC做分类) 对于SVM算法,我们首先 ...

季节性时间序列分析-SARIMAX模型的python实现