python時間序列分析包_時間序列分析之ARIMA上手-Python

airpassenger.csv：https://pan.baidu.com/s/1mi7dp6w

概念時間序列時間序列（或稱動态數列）是指将同一統計名額的數值按其發生的時間先後順序排列而成的數列。時間序列分析的主要目的是根據已有的曆史資料對未來進行預測。時間序列分析時間序列分析是根據系統觀察得到的時間序列資料，通過曲線拟合和參數估計來建立數學模型的理論和方法。時間序列分析常用于國民宏觀經濟控制、市場潛力預測、氣象預測、農作物害蟲災害預報等各個方面。

組成要素構成要素：長期趨勢，季節變動，循環變動，不規則變動.長期趨勢（ T ）現象在較長時期内受某種根本性因素作用而形成的總的變動趨勢

季節變動（ S ）現象在一年内随着季節的變化而發生的有規律的周期性變動

循環變動（ C ）現象以若幹年為周期所呈現出的波浪起伏形态的有規律的變動

不規則變動（I ）是一種無規律可循的變動，包括嚴格的随機變動和不規則的突發性影響很大的變動兩種類型

分析模型時間數列的組合模型加法模型：Y=T+S+C+I (Y,T 計量機關相同的總量名額)(S,C,I 對長期趨勢産生的或正或負的偏差)

乘法模型：Y=T·S·C·

4000

I(常用模型) (Y,T 計量機關相同的總量名額)(S,C,I 對原數列名額增加或減少的百分比)

ARIMA基本步驟擷取被觀測系統時間序列資料；

對資料繪圖，觀測是否為平穩時間序列；對于非平穩時間序列要先進行d階差分運算，化為平穩時間序列；

經過第二步處理，已經得到平穩時間序列。要對平穩時間序列分别求得其自相關系數ACF 和偏自相關系數PACF ，通過對自相關圖和偏自相關圖的分析，得到最佳的階層 p 和階數 q；

由以上得到的d、q、p ，得到ARIMA模型。然後開始對得到的模型進行模型檢驗。

Python ARIMA實踐資料：AirPassengers.csv使用的基礎庫： pandas,numpy,scipy,matplotlib,statsmodels。import pandas as pd

import numpy as np

import matplotlib.pylab as plt

from matplotlib.pylab import rcParams

from statsmodels.tsa.stattools import adfuller

from statsmodels.tsa.seasonal import seasonal_decompose

from statsmodels.tsa.stattools import acf, pacf

from statsmodels.tsa.arima_model import ARIMA

讀取資料：# data=pd.read_csv('/Users/wangtuntun/Desktop/AirPassengers.csv')

dateparse = lambda dates: pd.datetime.strptime(dates, '%Y-%m')

# paese_dates指定日期在哪列 ;index_dates将年月日的哪個作為索引 ;date_parser将字元串轉為日期

data = pd.read_csv('D:\\Competition\\AirPassengers.csv', parse_dates=['Month'], index_col='Month', date_parser=dateparse)

分析資料def test_stationarity(timeseries):

# 決定起伏統計

rolmean = pd.rolling_mean(timeseries, window=12) # 對size個資料進行移動平均

rol_weighted_mean = pd.ewma(timeseries, span=12) # 對size個資料進行權重移動平均

rolstd = pd.rolling_std(timeseries, window=12) # 偏離原始值多少

# 畫出起伏統計

orig = plt.plot(timeseries, color='blue', label='Original')

mean = plt.plot(rolmean, color='red', label='Rolling Mean')

weighted_mean = plt.plot(rol_weighted_mean, color='green', label='weighted Mean')

std = plt.plot(rolstd, color='black', label='Rolling Std')

plt.legend(loc='best')

plt.title('Rolling Mean & Standard Deviation')

plt.show(block=False)

# 進行df測試

print 'Result of Dickry-Fuller test'

dftest = adfuller(timeseries, autolag='AIC')

dfoutput = pd.Series(dftest[0:4], index=['Test Statistic', 'p-value', '#Lags Used', 'Number of observations Used'])

for key, value in dftest[4].items():

dfoutput['Critical value(%s)' % key] = value

print dfoutput

ts = data['#Passengers']

plt.plot(ts)

plt.show()

test_stationarity(ts)

plt.show()

可以看到長期趨勢和周期性的變動。# estimating

ts_log = np.log(ts)

# plt.plot(ts_log)

# plt.show()

moving_avg = pd.rolling_mean(ts_log, 12)

# plt.plot(moving_avg)

# plt.plot(moving_avg,color='red')

# plt.show()

ts_log_moving_avg_diff = ts_log - moving_avg

# print ts_log_moving_avg_diff.head(12)

ts_log_moving_avg_diff.dropna(inplace=True)

test_stationarity(ts_log_moving_avg_diff)

plt.show()

時間序列的差分d# 差分differencing

ts_log_diff = ts_log.diff(1)

ts_log_diff.dropna(inplace=True)

test_stationarity(ts_log_diff)

plt.show()

上圖看出一階差分大緻已經具有周期性，不妨繪制二階差分對比：ts_log_diff1 = ts_log.diff(1)

ts_log_diff2 = ts_log.diff(2)

ts_log_diff1.plot()

ts_log_diff2.plot()

plt.show()

基本已經沒有變化。是以使用一階差分。

分解decomposing# 分解decomposing

decomposition = seasonal_decompose(ts_log)

trend = decomposition.trend # 趨勢

seasonal = decomposition.seasonal # 季節性

residual = decomposition.resid # 剩餘的

plt.subplot(411)

plt.plot(ts_log,label='Original')

plt.legend(loc='best')

plt.subplot(412)

plt.plot(trend,label='Trend')

plt.legend(loc='best')

plt.subplot(413)

plt.plot(seasonal,label='Seasonarity')

plt.legend(loc='best')

plt.subplot(414)

plt.plot(residual,label='Residual')

plt.legend(loc='best')

plt.tight_layout()

plt.show()

預測

确定參數# 确定參數

lag_acf = acf(ts_log_diff, nlags=20)

lag_pacf = pacf(ts_log_diff, nlags=20, method='ols

c169

')

# q的擷取:ACF圖中曲線第一次穿過上置信區間.這裡q取2

plt.subplot(121)

plt.plot(lag_acf)

plt.axhline(y=0, linestyle='--', color='gray')

plt.axhline(y=-1.96 / np.sqrt(len(ts_log_diff)), linestyle='--', color='gray') # lowwer置信區間

plt.axhline(y=1.96 / np.sqrt(len(ts_log_diff)), linestyle='--', color='gray') # upper置信區間

plt.title('Autocorrelation Function')

# p的擷取:PACF圖中曲線第一次穿過上置信區間.這裡p取2

plt.subplot(122)

plt.plot(lag_pacf)

plt.axhline(y=0, linestyle='--', color='gray')

plt.axhline(y=-1.96 / np.sqrt(len(ts_log_diff)), linestyle='--', color='gray')

plt.axhline(y=1.96 / np.sqrt(len(ts_log_diff)), linestyle='--', color='gray')

plt.title('Partial Autocorrelation Function')

plt.tight_layout()

plt.show()

AR modelmodel = ARIMA(ts_log, order=(2, 1, 0))

result_AR = model.fit(disp=-1)

plt.plot(ts_log_diff)

plt.plot(result_AR.fittedvalues, color='red')

plt.title('AR model RSS:%.4f' % sum(result_AR.fittedvalues - ts_log_diff) ** 2)

plt.show()

MA model# MA model

model = ARIMA(ts_log, order=(0, 1, 2))

result_MA = model.fit(disp=-1)

plt.plot(ts_log_diff)

plt.plot(result_MA.fittedvalues, color='red')

plt.title('MA model RSS:%.4f' % sum(result_MA.fittedvalues - ts_log_diff) ** 2)

plt.show()

ARIMA model# ARIMA 将兩個結合起來 效果更好

model = ARIMA(ts_log, order=(2, 1, 2))

result_ARIMA = model.fit(disp=-1)

plt.plot(ts_log_diff)

plt.plot(result_ARIMA.fittedvalues, color='red')

plt.title('ARIMA RSS:%.4f' % sum(result_ARIMA.fittedvalues - ts_log_diff) ** 2)

plt.show()

predictions_ARIMApredictions_ARIMA_diff = pd.Series(result_ARIMA.fittedvalues, copy=True)

# print predictions_ARIMA_diff.head()#發現資料是沒有第一行的,因為有1的延遲

predictions_ARIMA_diff_cumsum = predictions_ARIMA_diff.cumsum()

# print predictions_ARIMA_diff_cumsum.head()

predictions_ARIMA_log = pd.Series(ts_log.ix[0], index=ts_log.index)

predictions_ARIMA_log = predictions_ARIMA_log.add(predictions_ARIMA_diff_cumsum, fill_value=0)

# print predictions_ARIMA_log.head()

predictions_ARIMA = np.exp(predictions_ARIMA_log)

plt.plot(ts)

plt.plot(predictions_ARIMA)

plt.title('predictions_ARIMA RMSE: %.4f' % np.sqrt(sum((predictions_ARIMA - ts) ** 2) / len(ts)))

plt.show()

文章出處：https://blog.csdn.net/shine19930820/article/details/72667656