时间序列预测模型实战案例(六)深入理解机器学习ARIMA包括差分和相关性分析_arima sklearn

​

import pandas as pd
import matplotlib.pyplot as plt
from pmdarima.arima import auto_arima
from math import sqrt
from sklearn.metrics import mean_squared_error
import statsmodels.tsa.api as smt
import statsmodels.api as sm
import numpy as np
 
firstDate = '2017-10-29'
 
def loadData(filterDate):
    """
    加载指定日期之后的数据，只加载第一天之后的数据，第一天可以指定
    :param filterDate:
    :return:
    """
    if filterDate < firstDate:
        filterDate = firstDate
    hp_all = pd.read_csv("hp.csv", index_col=['Date'], parse_dates=['Date'])
    hp_all = hp_all[filterDate:]
    return hp_all
 
 
def plotData(series):
    """
    对于输入数据进行折线图显示
    :param series:
    :return:
    """
    plt.rcParams['font.sans-serif'] = ['SimHei']
    plt.rcParams['axes.unicode_minus'] = False
    plt.figure(figsize=(15, 7))
    plt.plot(series)
    plt.title('原始数据')
    plt.grid(True)
    # plt.show()
 
 
def plotAllData(filterDate):
    """
    对指定日期之后的数据进行展示，初步处理是将数据中的空值替换为0
    :param filterDate:
    :return:
    """
    hp_all = loadData(filterDate)
    plotData(hp_all['Volume'].fillna(0, inplace=False))
 
 
def preProcessData(filterDate):
    """
    对指定日期之后的数据进行预处理，预处理方式将标记的数据使用周围数据的平均值替换
    :param filterDate:
    :return: hp 处理后的数据，其中标记的数据使用平滑值来代替；
    hp_rolling 窗口为17的平滑值，平滑方法是平均值；
    hp_raw 处理后的数据，标记的数据使用0来替换
    """
    hp_all = loadData(filterDate)
    hp = hp_all['Volume'].copy(deep=True)
    hp_all.loc[hp_all['IsNoted'] == 1, 'Volume'] = np.nan
    hp_raw = hp_all[['Volume']].fillna(0, inplace=False)
    hp_rolling = hp_all['Volume'].rolling(window=22, min_periods=1, center=True).mean()
    hp[hp_all['IsNoted'] == 1] = hp_rolling[hp_all['IsNoted'] == 1]
    return hp, hp_rolling, hp_raw
 
 
def plotProcessedData(filterDate):
    """
    对原始数据，平滑数据和处理后的数据进行展示
    :param filterDate:
    :return:
    """
    hp, hp_rolling, hp_raw = preProcessData(filterDate)
    plt.rcParams['font.sans-serif'] = ['SimHei']
    plt.rcParams['axes.unicode_minus'] = False
    plt.figure(figsize=(15, 7))
    plt.plot(hp_raw, label='Raw Data')
    # plt.plot(hp_rolling, label='Rolling Data')
    # plt.scatter(hp.index, hp, label='Processed Data', c='red', marker='v', s=5)
    plt.plot(hp, label='Processed Data')
    # plt.title('数据处理：原始数据，平滑数据，处理后的数据')
    plt.title('数据处理：原始数据，处理后的数据')
    plt.grid(True)
    plt.legend(loc='best')
    # plt.show()
 
 
def train(filterDate, div_para):
    hp, hp_rolling, hp_raw = preProcessData(filterDate)
    train_ds = hp[:int(div_para * len(hp))]
    test_ds = hp[int(div_para * len(hp)):]
 
    plt.rcParams['font.sans-serif'] = ['SimHei']
    plt.rcParams['axes.unicode_minus'] = False
    # plt.figure(figsize=(15, 7))
    # plt.title('训练数据和测试数据')
    # train_ds.plot()
    # test_ds.plot()
 
    model = auto_arima(train_ds, trace=True, error_action='ignore', suppress_warnings=True, seasonal=True,m=7,
                       stationary=False)
    model.fit(train_ds)
    forcast = model.predict(n_periods=len(test_ds))
    forcast = pd.DataFrame(forcast, index=test_ds.index, columns=['Prediction'])
    rms = sqrt(mean_squared_error(test_ds, forcast))
    plt.rcParams['font.sans-serif'] = ['SimHei']
    plt.rcParams['axes.unicode_minus'] = False
    plt.figure(figsize=(15, 7))
    # plt.figure('Forecast')
    plt.title('预测结果,参数：训练开始日期：{0}，共{1}天，所占百分比：{2}，评估均方根误差：{3}'.format(filterDate, len(train_ds), div_para, rms))
    plt.plot(train_ds, label='Train')
    plt.plot(test_ds, label='Test')
    plt.plot(forcast, label='Prediction')
    plt.legend()
    plt.show()
 
 
def tsplot(y, lags=None, figsize=(12, 7), style='bmh'):
    """
        Plot time series, its ACF and PACF, calculate Dickey–Fuller test
        y - timeseries
        lags - how many lags to include in ACF, PACF calculation
    """
    if not isinstance(y, pd.Series):
        y = pd.Series(y)
 
    with plt.style.context(style):
        fig = plt.figure(figsize=figsize)
        layout = (2, 2)
        ts_ax = plt.subplot2grid(layout, (0, 0), colspan=2)
        acf_ax = plt.subplot2grid(layout, (1, 0))
        pacf_ax = plt.subplot2grid(layout, (1, 1))
 
        y.plot(ax=ts_ax)
        p_value = sm.tsa.stattools.adfuller(y)[1]
        ts_ax.set_title('Time Series Analysis Plots\n Dickey-Fuller: p={0:.5f}'.format(p_value))
        smt.graphics.plot_acf(y, lags=lags, ax=acf_ax)
        smt.graphics.plot_pacf(y, lags=lags, ax=pacf_ax)
        plt.tight_layout()
 
 
 
 
if __name__ == '__main__':
    """
    "ARIMA",
    """
 
    startDate = '2015-01-15'
    # startDate = firstDate
    # hp_all = loadData(startDate)
    hp, hp_rolling, hp_raw = preProcessData(startDate)
    plotAllData(startDate)
    plotProcessedData(startDate)
 
    #  是否进行差分操作可视化
    tsplot_bool = False
    if tsplot_bool:
        ads = hp_raw
        # tsplot(ads.Ads, lags=60)
        ads_diff = ads.Volume - ads.Volume.shift(24)
        # tsplot(ads_diff[24:], lags=60)
        ads_diff = ads_diff - ads_diff.shift(1)
        tsplot(ads_diff[24 + 1:], lags=60)
 
    print("test")
    # 模型拟合以及预测
    startDate = firstDate
    train(startDate, 0.98)
 

def loadData(filterDate):
    """
    加载指定日期之后的数据，只加载第一天之后的数据，第一天可以指定
    :param filterDate:
    :return:
    """
    if filterDate < firstDate:
        filterDate = firstDate
    hp_all = pd.read_csv("hp.csv", index_col=['Date'], parse_dates=['Date'])
    hp_all = hp_all[filterDate:]
    return hp_all
 
 
def preProcessData(filterDate):
    """
    对指定日期之后的数据进行预处理，预处理方式将标记的数据使用周围数据的平均值替换
    :param filterDate:
    :return: hp 处理后的数据，其中标记的数据使用平滑值来代替；
    hp_rolling 窗口为17的平滑值，平滑方法是平均值；
    hp_raw 处理后的数据，标记的数据使用0来替换
    """
    hp_all = loadData(filterDate)
    hp = hp_all['Volume'].copy(deep=True)
    hp_all.loc[hp_all['IsNoted'] == 1, 'Volume'] = np.nan
    hp_raw = hp_all[['Volume']].fillna(0, inplace=False)
    hp_rolling = hp_all['Volume'].rolling(window=22, min_periods=1, center=True).mean()
    hp[hp_all['IsNoted'] == 1] = hp_rolling[hp_all['IsNoted'] == 1]
    return hp, hp_rolling, hp_raw

    model = auto_arima(train_ds, trace=True, error_action='ignore', suppress_warnings=True, seasonal=True, m=7,
                       stationary=False)
    model.fit(train_ds)

  forcast = model.predict(n_periods=len(test_ds))
    forcast = pd.DataFrame(forcast, index=test_ds.index, columns=['Prediction'])
    rms = sqrt(mean_squared_error(test_ds, forcast))

    plt.rcParams['font.sans-serif'] = ['SimHei']
    plt.rcParams['axes.unicode_minus'] = False
    plt.figure(figsize=(15, 7))
    # plt.figure('Forecast')
    plt.title('预测结果,参数：训练开始日期：{0}，共{1}天，所占百分比：{2}，评估均方根误差：{3}'.format(filterDate, len(train_ds), div_para, rms))
    plt.plot(train_ds, label='Train')
    plt.plot(test_ds, label='Test')
    plt.plot(forcast, label='Prediction')
    plt.legend()
    plt.show()