python数据分析之时间序列分析详情_python_

## 输出高清图像 %config InlineBackend.figure_format = 'retina' %matplotlib inline ## 图像显示中文的问题 import matplotlib matplotlib.rcParams['axes.unicode_minus']=False import seaborn as sns sns.set(font= "Kaiti",style="ticks",font_scale=1.4)

## 导入会使用到的相关库 import numpy as np import pandas as pd import matplotlib.pyplot as plt from statsmodels.tsa.stattools import * import statsmodels.api as sm import statsmodels.formula.api as smf from statsmodels.tsa.api import SimpleExpSmoothing,Holt,ExponentialSmoothing,AR,ARIMA,ARMA from statsmodels.graphics.tsaplots import plot_acf, plot_pacf import pmdarima as pm from sklearn.metrics import mean_absolute_error import pyflux as pf from fbprophet import Prophet ## 忽略提醒 import warnings warnings.filterwarnings("ignore")

## 读取时间序列数据,该数据包含:X1为飞机乘客数据,X2为一组随机数据 df = pd.read_csv("data/chap6/timeserise.csv") ## 查看数据的变化趋势 df.plot(kind = "line",figsize = (10,6)) plt.grid() plt.title("时序数据") plt.show()

## 白噪声检验Ljung-Box检验 ## 该检验用来检查序列是否为随机序列，如果是随机序列，那它们的值之间没有任何关系 ## 使用LB检验来检验序列是否为白噪声，原假设为在延迟期数内序列之间相互独立。 lags = [4,8,16,32] LB = sm.stats.diagnostic.acorr_ljungbox(df["X1"],lags = lags,return_df = True) print("序列X1的检验结果:\n",LB) LB = sm.stats.diagnostic.acorr_ljungbox(df["X2"],lags = lags,return_df = True) print("序列X2的检验结果:\n",LB) ## 如果P值小于0.05，说明序列之间不独立，不是白噪声 ''' 序列X1的检验结果: lb_stat lb_pvalue 4 427.738684 2.817731e-91 8 709.484498 6.496271e-148 16 1289.037076 1.137910e-264 32 1792.523003 0.000000e+00 序列X2的检验结果: lb_stat lb_pvalue 4 1.822771 0.768314 8 8.452830 0.390531 16 15.508599 0.487750 32 28.717743 0.633459 '''

## 序列的单位根检验，即检验序列的平稳性 dftest = adfuller(df["X2"],autolag='BIC') dfoutput = pd.Series(dftest[0:4], index=['adf','p-value','usedlag','Number of Observations Used']) print("X2单位根检验结果:\n",dfoutput) dftest = adfuller(df["X1"],autolag='BIC') dfoutput = pd.Series(dftest[0:4], index=['adf','p-value','usedlag','Number of Observations Used']) print("X1单位根检验结果:\n",dfoutput) ## 对X1进行一阶差分后的序列进行检验 X1diff = df["X1"].diff().dropna() dftest = adfuller(X1diff,autolag='BIC') dfoutput = pd.Series(dftest[0:4], index=['adf','p-value','usedlag','Number of Observations Used']) print("X1一阶差分单位根检验结果:\n",dfoutput) ## 一阶差分后 P值大于0.05, 小于0.1，可以认为其是平稳的 ''' X2单位根检验结果: adf -1.124298e+01 p-value 1.788000e-20 usedlag 0.000000e+00 Number of Observations Used 1.430000e+02 dtype: float64 X1单位根检验结果: adf 0.815369 p-value 0.991880 usedlag 13.000000 Number of Observations Used 130.000000 dtype: float64 X1一阶差分单位根检验结果: adf -2.829267 p-value 0.054213 usedlag 12.000000 Number of Observations Used 130.000000 dtype: float64 '''

## KPSS检验的原假设为：序列x是平稳的。 ## 对序列X2使用KPSS检验平稳性 dfkpss = kpss(df["X2"]) dfoutput = pd.Series(dfkpss[0:3], index=["kpss_stat"," p-value"," usedlag"]) print("X2 KPSS检验结果:\n",dfoutput) ## 接受序列平稳的原假设 ## 对序列X1使用KPSS检验平稳性 dfkpss = kpss(df["X1"]) dfoutput = pd.Series(dfkpss[0:3], index=["kpss_stat"," p-value"," usedlag"]) print("X1 KPSS检验结果:\n",dfoutput) ## 拒绝序列平稳的原假设 ## 对序列X1使用KPSS检验平稳性 dfkpss = kpss(X1diff) dfoutput = pd.Series(dfkpss[0:3], index=["kpss_stat"," p-value"," usedlag"]) print("X1一阶差分KPSS检验结果:\n",dfoutput) ## 接受序列平稳的原假设 ''' X2 KPSS检验结果: kpss_stat 0.087559 p-value 0.100000 usedlag 14.000000 dtype: float64 X1 KPSS检验结果: kpss_stat 1.052175 p-value 0.010000 usedlag 14.000000 dtype: float64 X1一阶差分KPSS检验结果: kpss_stat 0.05301 p-value 0.10000 usedlag 14.00000 dtype: float64 '''

## 检验ARIMA模型的参数d X1d = pm.arima.ndiffs(df["X1"], alpha=0.05, test="kpss", max_d=3) print("使用KPSS方法对序列X1的参数d取值进行预测,d = ",X1d) X1diffd = pm.arima.ndiffs(X1diff, alpha=0.05, test="kpss", max_d=3) print("使用KPSS方法对序列X1一阶差分后的参数d取值进行预测,d = ",X1diffd) X2d = pm.arima.ndiffs(df["X2"], alpha=0.05, test="kpss", max_d=3) print("使用KPSS方法对序列X2的参数d取值进行预测,d = ",X2d) ''' 使用KPSS方法对序列X1的参数d取值进行预测,d = 1 使用KPSS方法对序列X1一阶差分后的参数d取值进行预测,d = 0 使用KPSS方法对序列X1的参数d取值进行预测,d = 0 '''

## 检验SARIMA模型的参数季节阶数D X1d = pm.arima.nsdiffs(df["X1"], 12, max_D=2) print("对序列X1的季节阶数D取值进行预测,D = ",X1d) X1diffd = pm.arima.nsdiffs(X1diff, 12, max_D=2) print("序列X1一阶差分后的季节阶数D取值进行预测,D = ",X1diffd) ''' 对序列X1的季节阶数D取值进行预测,D = 1 序列X1一阶差分后的季节阶数D取值进行预测,D = 1 '''

## 对随机序列X2进行自相关和偏相关分析可视化 fig = plt.figure(figsize=(16,5)) plt.subplot(1,3,1) plt.plot(df["X2"],"r-") plt.grid() plt.title("X2序列波动") ax = fig.add_subplot(1,3,2) plot_acf(df["X2"], lags=60,ax = ax) plt.grid() ax = fig.add_subplot(1,3,3) plot_pacf(df["X2"], lags=60,ax = ax) plt.grid() plt.tight_layout() plt.show()

## 对非随机序列X1进行自相关和偏相关分析可视化 fig = plt.figure(figsize=(16,5)) plt.subplot(1,3,1) plt.plot(df["X1"],"r-") plt.grid() plt.title("X1序列波动") ax = fig.add_subplot(1,3,2) plot_acf(df["X1"], lags=60,ax = ax) plt.grid() ax = fig.add_subplot(1,3,3) plot_pacf(df["X1"], lags=60,ax = ax) plt.ylim([-1,1]) plt.grid() plt.tight_layout() plt.show()

## 对非随机序列X1一阶差分后的序列进行自相关和偏相关分析可视化 fig = plt.figure(figsize=(16,5)) plt.subplot(1,3,1) plt.plot(X1diff,"r-") plt.grid() plt.title("X1序列一阶差分后波动") ax = fig.add_subplot(1,3,2) plot_acf(X1diff, lags=60,ax = ax) plt.grid() ax = fig.add_subplot(1,3,3) plot_pacf(X1diff, lags=60,ax = ax) plt.grid() plt.tight_layout() plt.show()

## 时间序列的分解 ## 通过观察序列X1，可以发现其既有上升的趋势，也有周期性的趋势，所以可以将该序列进行分解 ## 使用乘法模型分解结果(通常适用于有增长趋势的序列) X1decomp = pm.arima.decompose(df["X1"].values,"multiplicative", m=12) ## 可视化出分解的结果 ax = pm.utils.decomposed_plot(X1decomp,figure_kwargs = {"figsize": (10, 6)}, show=False) ax[0].set_title("乘法模型分解结果") plt.show()

## 使用加法模型分解结果(通常适用于平稳趋势的序列) X1decomp = pm.arima.decompose(X1diff.values,"additive", m=12) ## 可视化出分解的结果 ax = pm.utils.decomposed_plot(X1decomp,figure_kwargs = {"figsize": (10, 6)}, show=False) ax[0].set_title("加法模型分解结果") plt.show()

## 数据准备 ## 对序列X1进行切分，后面的24个数据用于测试集 train = pd.DataFrame(df["X1"][0:120]) test = pd.DataFrame(df["X1"][120:]) ## 可视化切分后的数据 train["X1"].plot(figsize=(14,7), title= "乘客数量数据",label = "X1 train") test["X1"].plot(label = "X1 test") plt.legend() plt.grid() plt.show()

print(train.shape) print(test.shape) df["X1"].shape ''' (120, 1) (24, 1) (144,) '''

## 简单移动平均进行预测 y_hat_avg = test.copy(deep = False) y_hat_avg["moving_avg_forecast"] = train["X1"].rolling(24).mean().iloc[-1] ## 可视化出预测结果 plt.figure(figsize=(14,7)) train["X1"].plot(figsize=(14,7),label = "X1 train") test["X1"].plot(label = "X1 test"
                
提示：
                    本文由神整理自网络，如有侵权请联系本站删除！
                    

                    本站声明： 

                    1、本站所有资源均来源于互联网，不保证100%完整、不提供任何技术支持； 

                    2、本站所发布的文章以及附件仅限用于学习和研究目的;不得将用于商业或者非法用途；否则由此产生的法律后果，本站概不负责！
                
                

                

                                            
上一篇：Python Numpy中数组的集合操作详解_python_
                                                                
下一篇：python词云库wordcloud自定义词云制作步骤分享_python_