利用机器学习模型完成时间序列预测

1.简述

  • 时间序列数据是一种典型的数据,时间序列预测方法比较多。比如ARIMA模型、Prophet模型、指数平均法、滑动平均法等等。
  • 本文采用机器学习算法,如线性回归、随机森林等,完成时间序列预测,预测效果也比较好。

2.数据集

本文对应的数据集格式如下:

time value
2018-09-01 00:00 3221
2018-09-01 01:00 5515
2018-09-01 02:00 9971
.... ...
2018-09-05 01:00 4416
.... ...

如1小时粒度数据。根据历史数据,对未来一段时间数据进行预测。

3.特征介绍(精髓)

由于是单变量数据,特征主要包括两部分:平移特征和时间特征。
平移特征是指,将value向前已经平移操作shift。
时间特征指,每分钟均值、每小时均值、每工作日均值、是否节假日等。

平移特征
时间特征

当然,读者可以自行思考,是否有其他特征,也欢迎留言,大家一起探讨。

4.预测模型

import pandas as pd
import matplotlib.pylab as plt
import numpy as np
from sklearn.model_selection import cross_val_score
from sklearn.model_selection import TimeSeriesSplit
from sklearn.preprocessing import StandardScaler
from sklearn.linear_model import LassoCV, RidgeCV, LinearRegression
from sklearn.svm import SVR
from sklearn.ensemble import RandomForestRegressor
from sklearn.model_selection import GridSearchCV


# error计算误差
def mean_absolute_error(y_true, y_pred):
    return np.mean(np.abs((y_true - y_pred) / y_true)) * 100

# data split分割数据,训练集、测试集、预测集
def timeseries_split(x, y, test_size, pred_size):
    index_test = int(len(x) * (1 - test_size))
    x_train = x.iloc[:index_test]
    y_train = y.iloc[:index_test]
    x_test = x.iloc[index_test:len(x) - pred_size]
    y_test = y.iloc[index_test:len(x) - pred_size]
    x_pred = x.iloc[-pred_size:]
    y_pred = y.iloc[-pred_size:]
    return x_train, y_train, x_test, y_test, x_pred, y_pred

# calculate mean计算均值特征
def cal_mean(data, x_feature, y_feature):
    return dict(data.groupby(x_feature)[y_feature].mean())#利用分组,计算均值特征

# make feature计算平移特征
def build_feature(data, lag_start, lag_end, test_size, target_encoding=False, num_day_pred=1):#target_encoding是否要开启均值特征,num_day_pred预测多少天
    # build future data with 0
    last_date = data["time"].max()
    # 预测点个数,由数据粒度决定
    pred_points = int(num_day_pred * 24)  # 1h粒度,1day = 24个点
    pred_date = pd.date_range(start=last_date, periods=pred_points + 1, freq="1h")
    pred_date = pred_date[pred_date > last_date]  # 排除掉last_date(非预测), preiods = pred_points +1,也是因为last_date为非预测point,所以后延1个point
    future_data = pd.DataFrame({"time": pred_date, "y": np.zeros(len(pred_date))})#先将预测时间段的value设置为0,然后在利用shift和均值等构件特征,做预测
    # concat future data and last data
    df = pd.concat([data, future_data])
    df.set_index("time", drop=True, inplace=True)
    #print(df)
    # make feature
    # shift feature平移特征,lag_start,lag_end分别为shift平移多少,如从80-120,80,81,82,,,119,120.
    for i in range(lag_start, lag_end):
        df["lag_{}".format(i)] = df.y.shift(i)
    #diff feature#差分特征,对平移后的lag做差分操作,此特征作用不大
    df["diff_lag_{}".format(lag_start)] = df["lag_{}".format(lag_start)].diff(1)
    # time feature时间特征
    df["hour"] = df.index.hour
    # df["day"] = df.index.day
    # df["month"] = df.index.month
    df["minute"] = df.index.minute
    df["weekday"] = df.index.weekday
    df["weekend"] = df.weekday.isin([5, 6]) * 1
    df["holiday"] = 0
    df.loc["2018-10-01 00:00:00":"2018-10-07 23:00:00","holiday"] = 1
    #print(df)
    # df["holiday"]
    # average feature
    if target_encoding:  # 用test
        # 计算到已有数据截止,然后在映射到预测的数据中,这样就训练、测试、预测都有此特征
        df["weekday_avg"] = list(map(cal_mean(df[:last_date], "weekday", "y").get, df.weekday))#时间均值特征
        df["hour_avg"] = list(map(cal_mean(df[:last_date], "hour", "y").get, df.hour))
        df["weekend_avg"] = list(map(cal_mean(df[:last_date], "weekend", "y").get, df.weekend))
        df["minute_avg"] = list(map(cal_mean(df[:last_date], "minute", "y").get, df.minute))
        df = df.drop(["hour","minute","weekday", "weekend"], axis = 1)
    #one-hot没有作用
    #df = pd.get_dummies(df, columns = ["hour", "minute", "weekday", "weekend"])
    # data split
    y = df.dropna().y
    x = df.dropna().drop("y", axis=1)
    x_train, y_train, x_test, y_test, x_pred, y_pred = \
        timeseries_split(x, y, test_size=test_size, pred_size=pred_points)
    return x_train, y_train, x_test, y_test, x_pred, y_pred

# predict
def predict_future(model, scaler, x_pred, y_pred, lag_start, lag_end):#model拟合的模型,scaler归一化,lag平移特征,x_pred/y_pred预测x和y
    y_pred[0:lag_start] = model.predict(scaler.transform(x_pred[0:lag_start]))  # 预测到lag_start上一行
    for i in range(lag_start, len(x_pred)):
        last_line = x_pred.iloc[i-1]  # 已经预测数据的最后一行,还没预测数据的上一行,即shift,上一行填充到斜角下一行
        index = x_pred.index[i]
        x_pred.at[index, "lag_{}".format(lag_start)] = y_pred[i-1]
        x_pred.at[index, "diff_lag_{}".format(lag_start)] = y_pred[i-1] -  x_pred.at[x_pred.index[i-1], "lag_{}".format(lag_start)]
        for j in range(lag_start + 1, lag_end):  # 根据平移变换shift,前一个lag_{}列的值,shift后为下一个列的值
            x_pred.at[index, "lag_{}".format(j)] = last_line["lag_{}".format(j-1)]
        # 已经预测的最后一个值,赋值给lag_start对应的"lag_{}.format(lag_start)列
        # x_pred.at[index, "lag_{}".format(lag_start)] = y_pred[i - 1]
        y_pred[i] = model.predict(scaler.transform([x_pred.iloc[i]]))[0]
    return y_pred

# plot显示结果
def plot_result(y, y_fit, y_future):#y真实,y_fit拟合,y_future预测
    assert len(y) == len(y_fit)
    plt.figure(figsize=(16, 8))
    # plt.plot(y.index, y, "k.", label="y_orig")
    plt.plot(y.index, y, label="y_orig")
    plt.plot(y.index, y_fit, label="y_fit")
    plt.plot(y_future.index, y_future, "y", label="y_predict")
    error = mean_absolute_error(y, y_fit)
    plt.title("mean_absolute_error{0:.2f}%".format(error))
    plt.legend(loc="best")
    plt.grid(True)
    plt.show()

# coefs显示重要性
def plot_importance(model, x_train):
    coefs = pd.DataFrame(model.coef_, x_train.columns)
    #coefs = pd.DataFrame(model.feature_importances_, x_train.columns)
    coefs.columns = ["coefs"]
    coefs["coefs_abs"] = coefs.coefs.apply(np.abs)
    coefs = coefs.sort_values(by="coefs_abs", ascending=False).drop(["coefs_abs"], axis=1)
    plt.figure(figsize=(16, 6))
    coefs.coefs.plot(kind="bar")
    plt.grid(True, axis="y")
    plt.hlines(y=0, xmin=0, xmax=len(coefs), linestyles="dashed")
    plt.show()

# read data
if __name__ == "__main__":
    dataf = pd.read_csv("data.csv")
    dataf["time"] = pd.to_datetime(dataf["time"])
    dataf = dataf.sort_values("time")
    dataf.rename(columns={"sump": "y"}, inplace=True)
    lag_start = 80#要根据数据周期,调试
    lag_end = 120#平移特征
    x_train, y_train, x_test, y_test, x_pred, y_pred = build_feature \
        (dataf, lag_start=lag_start, lag_end=lag_end, test_size=0.3, target_encoding=True, num_day_pred=1)
    scaler = StandardScaler()
    x_train_scaled = scaler.fit_transform(x_train)
    x_test_scaled = scaler.transform(x_test)
    tscv = TimeSeriesSplit(n_splits=5)
    #lr = LassoCV(cv=tscv)
    lr = LinearRegression(normalize= "l1") 
    #可以尝试随机森林的效果,也不错。也可以做多模型结果融合,请自己尝试。
    #lr = RandomForestRegressor(n_estimators=100, max_depth=10) #lag_start = 288, lag_end = 320
    # lr = RidgeCV(cv = tscv)
    lr.fit(x_train_scaled, y_train)
    #train_score = lr.score(x_train_scaled, y_train)
    #test_score = lr.score(x_test_scaled, y_test)
    #print("num_tree", each, "score", train_score, test_score)
    # future预测
    y_future = predict_future(lr, scaler, x_pred, y_pred, lag_start, lag_end)
    #print(x_pred)
    # now 拟合
    y_fit = lr.predict(np.concatenate((x_train_scaled, x_test_scaled)))
    y = pd.concat([y_train, y_test])
    # 显示结果
    plt.figure(figsize=(16, 8))
    plt.plot(data["time"], data["sump"])
    plot_result(y, y_fit, y_future)
    y_future.to_csv("y_future_lr_test.csv")
    plot_importance(lr, x_train)

可以尝试随机森林等模型的效果。也可以做多模型结果融合,请自己尝试。

5.预测结果

预测结果

特征重要性如下:


特征重要性排序

介绍的比较简单,有不明白的,可以留言,探讨。

6.参考

主要参考了相关的博客,结合自己的数据,做了调整。
如果想学习,结合自己的业务,请见:
时间序列预测——深度好文中文
时间序列预测——ARIMA和随机森林对比英文
时间序列——各种方法分析英文
时间序列——预测方法对比中文
时间序列——Jason Brownlee博客系列推荐,有很多关于时间序列文章,但是英文的。

©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 194,088评论 5 459
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 81,715评论 2 371
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 141,361评论 0 319
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 52,099评论 1 263
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 60,987评论 4 355
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 46,063评论 1 272
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 36,486评论 3 381
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 35,175评论 0 253
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 39,440评论 1 290
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 34,518评论 2 309
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 36,305评论 1 326
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 32,190评论 3 312
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 37,550评论 3 298
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 28,880评论 0 17
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 30,152评论 1 250
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 41,451评论 2 341
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 40,637评论 2 335

推荐阅读更多精彩内容