利用ARIMA预测亚马逊商品排名

ARIMA算法简介：

ARIMA是一种基于统计学的时间序列算法。（Autoregressive Integrated Moving Average model），差分整合移动平均自回归模型，又称整合移动平均自回归模型（移动也可称作滑动），是时间序列预测分析方法之一。ARIMA(p，d，q)中，AR是“自回归”，p为自回归项数；MA为“滑动平均”，q为滑动平均项数，d为使之成为平稳序列所做的差分次数。

接下来代码算法评价指标使用的是：MAPE

MAPE 指平均绝对百分比误差，它是一种相对度量，它实际上将 MAD 尺度确定为百分比单位而不是变量的单位。平均绝对百分比误差是相对误差度量值，它使用绝对值来避免正误差和负误差相互抵消，您可以使用相对误差来比较各种时间序列模型预测的准确性。

 效果示意图：

 这里红线是真实的数据，我们选取某商品7天的排名变化

这里的蓝线是一步预测数据，非静态

这里的紫色线是动态预测的数据

我们可以看到一步预测相对于动态预测更为准确一些

 在上面的图中，我们可以观察到残差（真实值与预测值的差）是不相关的（右下图），并且没有表现出明显的周期性（左上图）。此外，残差和大致正态分布均为零（右上角图）。左下角的qq图显示，残差（蓝点）的有序分布与N（0，1）的标准正态分布样本的线性趋势大致一致。这再次有力地表明残差是正态分布的。

欢迎访问我的GitHub主页https://github.com/Huhaobangbangbang/time_series_project.git，上面有数据，以及随后将会放出LSTM，基于transformer的时间序列预测代码。

数据格式很简单，就是时间：排名

date,rank
"2022-10-11-11",6
"2022-10-11-12",3
"2022-10-11-13",4

代码如下：

# Import libraries
import warnings
import itertools
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import statsmodels.api as sm

# Defaults
plt.rcParams['figure.figsize'] = (20.0, 10.0)
plt.rcParams.update({'font.size': 12})
plt.style.use('ggplot')


# Load the data
data = pd.read_csv('/home/huhao/scripts/time_series/data/B06X6J3L65_rank.csv', engine='python', skipfooter=3)
# A bit of pre-processing to make it nicer
data['date']=pd.to_datetime(data['date'])
data.set_index(['date'], inplace=True)

# Plot the data
data.plot()
plt.ylabel('the product rank by date')
plt.xlabel('date')
plt.show()

# Define the d and q parameters to take any value between 0 and 1
q = d = range(0, 2)
# Define the p parameters to take any value between 0 and 3
p = range(0, 4)

# Generate all different combinations of p, q and q triplets
pdq = list(itertools.product(p, d, q))

# Generate all different combinations of seasonal p, q and q triplets
seasonal_pdq = [(x[0], x[1], x[2], 12) for x in list(itertools.product(p, d, q))]

print('Examples of parameter combinations for Seasonal ARIMA...')
print('SARIMAX: {} x {}'.format(pdq[1], seasonal_pdq[1]))
print('SARIMAX: {} x {}'.format(pdq[1], seasonal_pdq[2]))
print('SARIMAX: {} x {}'.format(pdq[2], seasonal_pdq[3]))
print('SARIMAX: {} x {}'.format(pdq[2], seasonal_pdq[4]))


train_data = data['2022-10-11-11':'2022-10-22-23']
test_data = data['2022-10-23-00':'2022-10-24-10']

warnings.filterwarnings("ignore") # specify to ignore warning messages

AIC = []
SARIMAX_model = []
for param in pdq:
    for param_seasonal in seasonal_pdq:
        try:
            mod = sm.tsa.statespace.SARIMAX(train_data,
                                            order=param,
                                            seasonal_order=param_seasonal,
                                            enforce_stationarity=False,
                                            enforce_invertibility=False)

            results = mod.fit()

            #print('SARIMAX{}x{} - AIC:{}'.format(param, param_seasonal, results.aic), end='\r')
            AIC.append(results.aic)
            SARIMAX_model.append([param, param_seasonal])
        except:
            continue  


print('The smallest AIC is {} for model SARIMAX{}x{}'.format(min(AIC), SARIMAX_model[AIC.index(min(AIC))][0],SARIMAX_model[AIC.index(min(AIC))][1]))
# Let's fit this model
mod = sm.tsa.statespace.SARIMAX(train_data,
                                order=SARIMAX_model[AIC.index(min(AIC))][0],
                                seasonal_order=SARIMAX_model[AIC.index(min(AIC))][1],
                                enforce_stationarity=False,
                                enforce_invertibility=False)

results = mod.fit()


results.plot_diagnostics(figsize=(20, 14))
plt.show()

pred0 = results.get_prediction(start='2022-10-21-00', dynamic=False)
pred0_ci = pred0.conf_int()
pred1 = results.get_prediction(start='2022-10-21-00', dynamic=True)
pred1_ci = pred1.conf_int()
pred2 = results.get_forecast(31)
pred2_ci = pred2.conf_int()
print(pred2_ci)
ax = data.plot(figsize=(20, 16))
pred0.predicted_mean.plot(ax=ax, label='1-step-ahead Forecast (get_predictions, dynamic=False)')
pred1.predicted_mean.plot(ax=ax, label='Dynamic Forecast (get_predictions, dynamic=True)')
# pred2.predicted_mean.plot(ax=ax, label='Dynamic Forecast (get_forecast)')
# ax.fill_between(pred2_ci.index, pred2_ci.iloc[:, 0], pred2_ci.iloc[:, 1], color='k', alpha=.1)
plt.ylabel('Monthly airline passengers (x1000)')
plt.xlabel('Date')
plt.legend()
plt.show()
prediction = pred2.predicted_mean.values
# flatten nested list
truth = list(itertools.chain.from_iterable(test_data.values))
# Mean Absolute Percentage Error
MAPE = np.mean(np.abs((truth - prediction) / truth)) * 100

print('The Mean Absolute Percentage Error for the forecast of product rank is {:.2f}%'.format(MAPE))

# 参考链接 https://github.com/gmonaci/ARIMA.git