周期性时间序列的预测
背景
最近在研究时间序列的时候，发现很多序列具有很强的周期性，那如何对此类序列进行预测呢？
数据处理
挑选一个如下图的具有周期性的时间序列。该序列是取得是过去7天的数据，每小时一个点，一共7*24个点。
划分数据集
我们取前六天的数据做训练，第七天做测试集。
平滑处理
时间序列经常会出现毛刺的点，需要做平滑处理才能分析，类似上图中的数据。消除数据的毛刺，可以用移动平均法，但是移动平均有时候处理完后并不能使数据平滑，我这里采用的方法很简单，但效果还不错：把每个点与上一点的变化值作为一个新的序列，对这里边的异常值，也就是变化比较离谱的值剃掉，用前后数据的均值填充：
def diff_smooth(ts, interval):
'''时间序列平滑处理'''
# 间隔为1小时
wide = interval/60
# 差分序列
dif = ts.diff().dropna()
# 描述性统计得到：min，25%，50%，75%，max值
td = dif.describe()
# 定义高点阈值，1.5倍四分位距之外
high = td['75%'] + 1.5 * (td['75%'] - td['25%'])
# 定义低点阈值
low = td['25%'] - 1.5 * (td['75%'] - td['25%'])
i = 0
forbid_index = dif[(dif > high) | (dif < low)].index
while i < len(forbid_index) - 1:
# 发现连续多少个点变化幅度过大
n = 1
# 异常点的起始索引
start = forbid_index[i]
while forbid_index[i+n] == start + datetime.timedelta(minutes=n):
n += 1
i += n - 1
# 异常点的结束索引
end = forbid_index[i]
# 用前后值的中间值均匀填充
value = np.linspace(ts[start - datetime.timedelta(minutes=wide)], ts[end + datetime.timedelta(minutes=wide)], n)
ts[start: end] = value
i += 1
return ts
经过处理以后，上图的时间序列得到了平滑处理，效果如下图。
周期性分解
具有周期性特征的序列需要将周期性特征提取出来。python里面的statsmodels工具包里面有针对周期性分解的函数seasonal_decompose，我们可以将序列进行分解。seasonal_decompose这个函数里面有个two_sided的参数，默认是True。Trend处理的时候用到移动平均的方法，熟悉此方法的读者就会发现，经过该方法处理以后，序列收尾两段有一部分数据缺失了，但是如果该参数为FALSE，则只有开始的时候有一段缺失值。
decomposition = seasonal_decompose(smooth_data, two_sided=False)
trend = decomposition.trend
seasonal = decomposition.seasonal
residual = decomposition.resid
图3中的第一张图是observed，体现的原始数据；第二张是trend，体现的是分解出来的趋势部分；第三张是seasonal，体现的是周期部分；最后是residual，体现的是残差部分。
本文采用的是seasonal_decompose的加法模型进行的分解，即 observed = trend + seasonal + residual，另还有乘法模型。在建模的时候，只针对trend部分学习和预测，如何将trend的预测结果加工成合理的最终结果？后面会有介绍。
预测
我们对trend部分进行预测，最后再加上seasonal部分。对trend的预测，我们采用ARIMA模型。熟悉该模型的都知道，需要确定三个参数p，q和d，可以使用aic和bic的方法进行定阶，可以查阅相关的文献。
trend.dropna(inplace=True)
# arima的训练参数order =（p,d,q），具体意义查看官方文档，调参过程略。
predict_model = ARIMA(trend, (1,1,2), freq='H').fit(disp=-1, method='css')
得到模型以后，就可以进行预测。
# 续接train，生成长度为n的时间索引，赋给预测序列
n = 24
predict_time_index = pd.date_range(start=trend.index[-1], periods=n+1, freq='H')[1:]
trend_predict = predict_model.forecast(n)[0]
'''
为预测出的趋势数据添加周期数据和残差数据
'''
values = []
low_conf_values = []
high_conf_values = []
# enumerate() 函数用于将一个可遍历的数据对象(如列表、元组或字符串)组合为一个索引序列,同时列出数据和数据下标,一般用在for循环当中。
for i, t in enumerate(predict_time_index):
trend_part = trend_predict[i]
# 相同时间点的周期数据均值
# t为2018-08-09 15:18:00类型的时间,t.time()为15:18:00类型的时间
season_part = seasonal[seasonal.index.time == t.time()].mean()
# 趋势 + 周期 + 误差界限
predict = trend_part + season_part
values.append(predict)
# 得到预测值
final_predict = pd.Series(values, index=predict_time_index, name='predict')
下面是预测的结果，从图中可以看到预测的结果将周期性的特征完美地体现出来了。
评估
对第七天作出预测，评估的指标为均方根误差rmse，本序列的rmse小于5，效果还是不错的。
总结
本文介绍了周期性序列的预测方法，你可能会问并不是所有的序列都具有周期性，事实确实如此，接下来几篇博客，我会重点介绍周期性检测的一些方法。希望此博客对您研究时间序列有所帮助。
