【销售预测 ARIMA模型】ARIMA模型预测每天的销售额

本文主要是介绍【销售预测 ARIMA模型】ARIMA模型预测每天的销售额，希望对大家解决编程问题提供一定的参考价值，需要的开发者们随着小编来一起学习吧！
输入数据txt格式：
2017-05-01 100
2017-05-02 200
…….
python 实现arima：
# encoding: utf-8"""
function：时间序列预测ARIMA模型预测每天的销售额
author:dongli
date:2018-05-25
"""# 导入库
import numpy as np  # numpy库
from statsmodels.graphics.tsaplots import plot_acf, plot_pacf  # acf和pacf展示库
from statsmodels.tsa.stattools import adfuller  # adf检验库
from statsmodels.stats.diagnostic import acorr_ljungbox  # 随机性检验库
from statsmodels.tsa.arima_model import ARMA  # ARMA库
import matplotlib.pyplot as plt  # matplotlib图形展示库
import prettytable  # 导入表格库
import pandas as pd  #导入数据处理库
from datetime import datetime,timedelta #导入时间库
import warnings
import time
time1=time.time()
from sklearn.externals import joblib  #模型导出
# 忽略提醒
warnings.filterwarnings("ignore")# 把字符串转成datetime
def string_toDatetime(string):''':param string: 日期字符串:return: 返回日期格式数据'''return datetime.strptime(string, "%Y-%m-%d")# 定义一个日期列表
def datelist(beginDate, endDate):''':beginDate string: 日期字符串:endDate string: 日期字符串:return: 返回日期格式列表'''# beginDate, endDate是形如‘20160601’的字符串或datetime 2018-04-23格式date_l=[datetime.strftime(x,'%Y-%m-%d') for x in list(pd.date_range(start=beginDate, end=endDate))]return date_l# 多次用到的表格
def pre_table(table_name, table_rows):''':param table_name: 表格名称，字符串列表:param table_rows: 表格内容，嵌套列表:return: 展示表格对象'''table = prettytable.PrettyTable()  # 创建表格实例table.field_names = table_name  # 定义表格列名for i in table_rows:  # 循环读多条数据table.add_row(i)  # 增加数据return table# 数据平稳处理
def get_best_log(ts, max_log=5, rule1=True, rule2=True):''':param ts: 时间序列数据，Series类型:param max_log: 最大log处理的次数，int型:param rule1: rule1规则布尔值，布尔型:param rule2: rule2规则布尔值，布尔型:return: 达到平稳处理的最佳次数值和处理后的时间序列'''if rule1 and rule2:  # 如果两个规则同时满足return 0, ts  # 直接返回0和原始时间序列数据else:  # 只要有一个规则不满足for i in range(1, max_log):  # 循环做log处理ts = np.log(ts)  # log处理adf, pvalue1, usedlag, nobs, critical_values, icbest = adfuller(ts)  # 稳定性（ADF）检验lbvalue, pvalue2 = acorr_ljungbox(ts, lags=1)  # 白噪声（随机性）检验rule_1 = (adf < critical_values['1%'] and adf < critical_values['5%'] and adf < critical_values['10%'] and pvalue1 < 0.01)  # 稳定性（ADF）检验规则rule_2 = (pvalue2 < 0.05)  # 白噪声（随机性）规则rule_3 = (i < 5)if rule_1 and rule_2 and rule_3:  # 如果同时满足条件# print ('The best log n is: {0}'.format(i))  # 打印输出最佳次数return i, ts  # 返回最佳次数和处理后的时间序列# 还原经过平稳处理的数据
def recover_log(ts, log_n):''':param ts: 经过log方法平稳处理的时间序列，Series类型:param log_n: log方法处理的次数，int型:return: 还原后的时间序列'''for i in range(1, log_n + 1):  # 循环多次ts = np.exp(ts)  # log方法还原return ts  # 返回时间序列# 稳定性（ADF）检验
def adf_val(ts, ts_title, acf_title, pacf_title):''':param ts: 时间序列数据，Series类型:param ts_title: 时间序列图的标题名称，字符串:param acf_title: acf图的标题名称，字符串:param pacf_title: pacf图的标题名称，字符串:return: adf值、adf的p值、三种状态的检验值'''# plt.figure()# plt.plot(ts)  # 时间序列图# plt.title(ts_title)  # 时间序列标题# plt.show()# plot_acf(ts, lags=20, title=acf_title).show()  # 自相关检测# plot_pacf(ts, lags=20, title=pacf_title).show()  # 偏相关检测adf, pvalue, usedlag, nobs, critical_values, icbest = adfuller(ts)  # 稳定性（ADF）检验table_name = ['adf', 'pvalue', 'usedlag', 'nobs', 'critical_values', 'icbest']  # 表格列名列表table_rows = [[adf, pvalue, usedlag, nobs, critical_values, icbest]]  # 表格行数据，嵌套列表adf_table = pre_table(table_name, table_rows)  # 获得平稳性展示表格对象# print ('stochastic score')  # 打印标题# print (adf_table)  # 打印展示表格return adf, pvalue, critical_values,  # 返回adf值、adf的p值、三种状态的检验值# 白噪声（随机性）检验
def acorr_val(ts):''':param ts: 时间序列数据，Series类型:return: 白噪声检验的P值和展示数据表格对象'''lbvalue, pvalue = acorr_ljungbox(ts, lags=1)  # 白噪声检验结果table_name = ['lbvalue', 'pvalue']  # 表格列名列表table_rows = [[lbvalue, pvalue]]  # 表格行数据，嵌套列表acorr_ljungbox_table = pre_table(table_name, table_rows)  # 获得白噪声检验展示表格对象# print ('stationarity score')  # 打印标题# print (acorr_ljungbox_table)  # 打印展示表格return pvalue  # 返回白噪声检验的P值和展示数据表格对象# arma最优模型训练
def arma_fit(ts):''':param ts: 时间序列数据，Series类型:return: 最优状态下的p值、q值、arma模型对象、pdq数据框和展示参数表格对象'''max_count = int(len(ts) / 10)  # 最大循环次数最大定义为记录数的10%bic = float('inf')  # 初始值为正无穷tmp_score = []  # 临时p、q、aic、bic和hqic的值的列表for tmp_p in range(max_count + 1):  # p循环max_count+1次for tmp_q in range(max_count + 1):  # q循环max_count+1次model = ARMA(ts, order=(tmp_p, tmp_q))  # 创建ARMA模型对象try:results_ARMA = model.fit(disp=-1, method='css')  # ARMA模型训练except:continue  # 遇到报错继续finally:tmp_aic = results_ARMA.aic  # 模型的获得aictmp_bic = results_ARMA.bic  # 模型的获得bictmp_hqic = results_ARMA.hqic  # 模型的获得hqictmp_score.append([tmp_p, tmp_q, tmp_aic, tmp_bic, tmp_hqic])  # 追加每个模型的训练参数和结果if tmp_bic < bic:  # 如果模型bic小于最小值，那么获得最优模型ARMA的下列参数：p = tmp_p  # 最优模型ARMA的p值q = tmp_q  # 最优模型ARMA的q值model_arma = results_ARMA  # 最优模型ARMA的模型对象aic = tmp_bic  # 最优模型ARMA的aicbic = tmp_bic  # 最优模型ARMA的bichqic = tmp_bic  # 最优模型ARMA的hqicpdq_metrix = np.array(tmp_score)  # 将嵌套列表转换为矩阵pdq_pd = pd.DataFrame(pdq_metrix, columns=['p', 'q', 'aic', 'bic', 'hqic'])  # 基于矩阵创建数据框table_name = ['p', 'q', 'aic', 'bic', 'hqic']  # 表格列名列表table_rows = [[p, q, aic, bic, hqic]]  # 表格行数据，嵌套列表parameter_table = pre_table(table_name, table_rows)  # 获得最佳ARMA模型结果展示表格对象# print ('each p/q traning record')  # 打印标题# print (pdq_pd)  # 打印输出每次ARMA拟合结果，包含p、d、q以及对应的AIC、BIC、HQIC# print ('best p and q')  # 打印标题# print (parameter_table)  # 输出最佳ARMA模型结果展示表格对象return model_arma  # 最优状态下的arma模型对象# 模型训练和效果评估
def train_test(model_arma, ts, log_n, rule1=True, rule2=True):''':param model_arma: 最优ARMA模型对象:param ts: 时间序列数据，Series类型:param log_n: 平稳性处理的log的次数，int型:param rule1: rule1规则布尔值，布尔型:param rule2: rule2规则布尔值，布尔型:return: 还原后的时间序列'''train_predict = model_arma.predict()  # 得到训练集的预测时间序列if not (rule1 and rule2):  # 如果两个条件有任意一个不满足train_predict = recover_log(train_predict, log_n)  # 恢复平稳性处理前的真实时间序列值ts = recover_log(ts, log_n)  # 时间序列还原处理ts_data_new = ts[train_predict.index]  # 将原始时间序列数据的长度与预测的周期对齐RMSE = np.sqrt(np.sum((train_predict - ts_data_new) ** 2) / ts_data_new.size)  # 求RMSEprint("均方根误差为:%s" %RMSE)# # 对比训练集的预测和真实数据# plt.figure()  # 创建画布# train_predict.plot(label='predicted data', style='--')  # 以虚线展示预测数据# ts_data_new.plot(label='raw data')  # 以实线展示原始数据# plt.legend(loc='best')  # 设置图例位置# plt.title('raw data and predicted data with RMSE of %.2f' % RMSE)  # 设置标题# plt.show()  # 展示图像return ts  # 返回还原后的时间序列# 预测未来指定时间项的数据
def predict_data(model_arma, ts, log_n,start, end,rule1=True, rule2=True):''':param model_arma: 最优ARMA模型对象:param ts: 时间序列数据，Series类型:param log_n: 平稳性处理的log的次数，int型:param start: 要预测数据的开始时间索引:param end: 要预测数据的结束时间索引:param rule1: rule1规则布尔值，布尔型:param rule2: rule2规则布尔值，布尔型:return: 无'''m1=string_toDatetime(end)-string_toDatetime(start)# 预测未来指定时间项的数据predict_ts = model_arma.predict(start=len(ts)-5, end=len(ts)+m1.days+1)beginDate =ts.index[len(ts)-6]endDate = string_toDatetime(end)m = datelist(beginDate, endDate)# print(m)print ('-----------predict data----------')  # 打印标题if not (rule1 and rule2):  # 如果两个条件有任意一个不满足predict_ts = recover_log(predict_ts, log_n)  # 还原数据predict_ts_dataframe=pd.DataFrame(predict_ts,columns = ['number'])# print(predict_ts_dataframe)predict_ts_dataframe['day']=mpre_result=predict_ts_dataframe[predict_ts_dataframe['day']>=start]kk=list(pre_result.loc[:,'number'].apply(lambda x: x*10))# 展示预测数据pre_result2=pd.DataFrame({"day":pre_result['day'],"number":kk})print(pre_result2)#####输出csvpre_result2.to_csv("C:\\Users\\xiaohu\\Desktop\\销售数据预测\\训练数据\\pre_result2.csv", index=False)if __name__ == '__main__':# 读取数据# date_parse = lambda dates: pd.datetime.strptime(dates, '%Y/%m/%d')  # 创建解析列的功能对象# df = pd.read_table('C:\\Users\\xiaohu\\Desktop\销售数据预测\\训练数据\\train_2018.txt', delimiter='\t', index_col='date', date_parser=date_parse)  # 读取数据# 读取数据date_parse = lambda dates: pd.datetime.strptime(dates, '%Y-%m-%d')  # 创建解析列的功能对象df = pd.read_table('C:\\Users\\xiaohu\\Desktop\\book\\chapter4\\time_series4.txt', delimiter='\t', index_col='date',date_parser=date_parse)  # 读取数据# print(df)ts_data = df['number'].astype('float32')  # 将列转换为float32类型# print ('data summary')  # 打印标题# print (ts_data.describe())  # 打印输出时间序列数据概况# 原始数据检验adf, pvalue1, critical_values = adf_val(ts_data, 'raw time series', 'raw acf', 'raw pacf')  # 稳定性检验pvalue2 = acorr_val(ts_data)  # 白噪声检验# 创建用于区分是否进行平稳性处理的规则rule1 = (adf < critical_values['1%'] and adf < critical_values['5%'] and adf < critical_values['10%'] and pvalue1 < 0.01)  # 稳定性检验rule2 = (pvalue2[0,] < 0.05)  # 白噪声检验# 对时间序列做稳定性处理log_n, ts_data = get_best_log(ts_data, max_log=5, rule1=rule1, rule2=rule2)# 再次做检验adf, pvalue1, critical_values = adf_val(ts_data, 'final time series', 'final acf', 'final pacf')  # 稳定性检验pvalue2 = acorr_val(ts_data)  # 白噪声检验# 训练最佳ARMA模型并输出相关参数和对象# model_arma = arma_fit(ts_data)# sklearn中提供了高效的模型持久化模块joblib，将模型保存至硬盘。#joblib.dump(model_arma, 'C:\\Users\\xiaohu\\Desktop\\销售数据预测\\训练数据\\model_arma_day.model')# 模型导出model_arma = joblib.load('C:\\Users\\xiaohu\\Desktop\\销售数据预测\\训练数据\\model_arma_day.model')# 模型导入## 模型训练和效果评估ts_data = train_test(model_arma, ts_data, log_n, rule1=rule1, rule2=rule2)# print(ts_data)# # 模型预测应用start ='2018-06-04'  # 设置预测开始的时间索引end ='2018-06-10'  # 设置预测结束的时间索引predict_data(model_arma, ts_data, log_n,start,end ,rule1=rule1, rule2=rule2)  # 预测时间序列数据time2=time.time()print('总共耗时：' + str(time2 - time1) + 's')