Energy Consumption (能源消耗）数据集预测（1）：1小时从听说prophet 到完成预测（含安装）

数据如琥珀
1,727 阅读8分钟

前言

上个专辑我们分析了COVID19的数据，主要从可视化的角度进行了分析。 其实 kaggle上COVID19的数据集多达上千个(链接)。 有兴趣的可以下载更多的内容研究，后续我也会考虑使用其中一部分数据进行建模。

今天我们来入手另一个经典数据集，Energy Consumption (能源消耗）。该数据集是美国PJM公司（相当于国家电网或者南方电网）统计的能源消耗情况，时间跨度为达20年，数据集的采样频率为小时。

时间序列预测简述

时间序列，一般指一维的时间序列（更准确的叫一元时间序列），这是约定俗成的省略。由于也有多元时间序列，所以在和别人交流时，需要判断ta所说的时间序列是不是特指一元时间序列。

我遇到的时间序列：

  • 工业控制中的几乎所有信号（measurement），都是一个随时间变化的数据。
  • 股市的股票
  • 几乎人类所有的可量化活动，都可以画在时间轴上

时间序列的预测就是用过去的信息，预测未来的信息。

  • 股票明天涨还是涨呢？（调皮）
  • 如何根据特朗普的民调支持率预测下个月的支持率是多少呢？
  • 如何预测明天的天气温度呢？

prophet 简述

你可能听说过prophet，隐约有如下印象：

  • facebook公司近几年开开源的算法
  • 专门针对时间序列预测
  • 听说能考虑节假日等的影响，很牛 。。。

够了，知道这些就可以阅读本文了。至于它如何分解时间序列，如何划分changepoint，等等。实战过程中你自然会想着去了解。为啥呢？想提升模型的效果，总得尝试调参啊。今天我就花了一个小时，完成从安装prophet（包含入坑出坑）到模型预测。

prophet 安装

这里特地花一点时间来说一下prophet的安装过程，因为它里面有坑，坑里面还有天坑。否则就是半小时从听说prophet 到完成预测了。

先说一下本人这次使用的环境：CPU only的笔记本+ pycharm 软件。

第一坑

打开prophet的官网，在安装指导页面，你会被这简单而又熟悉的pip代码吸引。

结果，安装失败。

回去仔细阅读安装说明，下面还有小字：

The major dependency that Prophet has is pystan. PyStan has its own installation instructions. Install pystan with pip before using pip to install fbprophet.

第二坑

pystan 一眼就看出是可以pip安装的。于是pip install pystan，没有报错，心里窃喜。然后pip install fbprophet

结果。。。还是失败！！

回去仔细阅读pystan的说明：

天坑

原来pystan也有一堆东西要提前安装，但是我一看到要安装c++，我瑟瑟发抖。因为这中间的坑不是一篇文章能说完的。再者，我不想弄乱自己windows的系统。

避免坑

于是花了一点时间仔细看说明，想要节省时间，就用pystan安装手册的推荐的用Anaconda 安装吧。 于是我的安装步骤如下：

  • 安装minicoda（洁癖），这就是典型的windows 软件安装方式：运行-->下一步-->完成
  • 在pycharm中配置采用anaconda 作为python的环境
  • conda install anaconda
  • conda install pystan
  • conda install -c conda-forge fbprophet

终于成功了，赶紧码代码。 一如既往，代码优先。

准备数据

从kaggle上下载数据，放到本地文件夹data中。 首先导入模块，包括数据分析经典的pandas，sklearn，numpy，matplotlib 以及我们本文的明星Prophet。

from fbprophet import Prophet
import pandas as pd
import matplotlib.pyplot as plt
from sklearn.metrics import mean_absolute_error
import numpy as np

简单说一下数据集，kaggle提供了PJM 各个区域（客户群）的数据。我们快速浏览每个数据的columns，发现除了pjm_hourly_est文件，剩下的都是两列数据，一个时间列，一个耗电量（MW）。pjm_hourly_est 不过是把所有其他文件的信息合并到一个文件而已。

AEP_hourly
Index(['Datetime', 'AEP_MW'], dtype='object')
******************************
COMED_hourly
Index(['Datetime', 'COMED_MW'], dtype='object')
******************************
DAYTON_hourly
Index(['Datetime', 'DAYTON_MW'], dtype='object')
******************************
DEOK_hourly
Index(['Datetime', 'DEOK_MW'], dtype='object')
******************************
DOM_hourly
Index(['Datetime', 'DOM_MW'], dtype='object')
******************************
DUQ_hourly
Index(['Datetime', 'DUQ_MW'], dtype='object')
******************************
EKPC_hourly
Index(['Datetime', 'EKPC_MW'], dtype='object')
******************************
FE_hourly
Index(['Datetime', 'FE_MW'], dtype='object')
******************************
NI_hourly
Index(['Datetime', 'NI_MW'], dtype='object')
******************************
PJME_hourly
Index(['Datetime', 'PJME_MW'], dtype='object')
******************************
PJMW_hourly
Index(['Datetime', 'PJMW_MW'], dtype='object')
******************************
pjm_hourly_est
Index(['Datetime', 'AEP', 'COMED', 'DAYTON', 'DEOK', 'DOM', 'DUQ', 'EKPC',
       'FE', 'NI', 'PJME', 'PJMW', 'PJM_Load'],
      dtype='object')
******************************
PJM_Load_hourly
Index(['Datetime', 'PJM_Load_MW'], dtype='object')
******************************

我们就是想看看prophet如何预测，所以就使用第一个数据集吧。

file = 'data/AEP_hourly.csv'
raw_df = pd.read_csv(file)
print(raw_df.head(5))
print(raw_df.tail(5))

结果如下，看来时间跨度从2004年末到2018年初。

              Datetime   AEP_MW
0  2004-12-31 01:00:00  13478.0
1  2004-12-31 02:00:00  12865.0
2  2004-12-31 03:00:00  12577.0
3  2004-12-31 04:00:00  12517.0
4  2004-12-31 05:00:00  12670.0
                   Datetime   AEP_MW
121268  2018-01-01 20:00:00  21089.0
121269  2018-01-01 21:00:00  20999.0
121270  2018-01-01 22:00:00  20820.0
121271  2018-01-01 23:00:00  20415.0
121272  2018-01-02 00:00:00  19993.0

调用prophet建模

对于prophet，啥都不知道，那就看指导手册吧，每个指导手册都会有quick start，让你快速调用。

数据稍作处理

对于prophet，输入的数据需要有以下格式：

  • Dataframe (我们已经满足)
  • 两列，第一列时间，第二列是观察对象（我们满足）
  • 第一列数据格式为YYYY-MM-DD HH:MM:SS 的Timestamp (我们不满足，因为pandas 读到的是默认字符串)
  • 第一列名称要为'ds',第二列名称为'y' (我们不满足)

所以，根据上述要求，我们来处理数据：

dt_format = '%Y-%m-%d %H:%M:%S'  # 处理时间格式
raw_df['Datetime']= pd.to_datetime(raw_df['Datetime'],format=dt_format)
raw_df.rename({'AEP_MW':'y','Datetime':'ds'},axis=1,inplace=True) # 重命名

print(raw_df.info())

Data columns (total 2 columns):
 #   Column  Non-Null Count   Dtype         
---  ------  --------------   -----         
 0   ds      121273 non-null  datetime64[ns]
 1   y       121273 non-null  float64       
dtypes: datetime64[ns](1), float64(1)

训练和测试数据集划分

对于机器学习，划分训练和测试数据集是必须的。因为时间序列有时间上的连贯性，不能轻易的用sklearn的train_test_split 函数划分。 这里我们划分很简单，我们刚才已经预览了数据集的首尾，知道器时间跨度为2004 到2018，有十几年的数据，所以我们这里就划分为最后一年为测试集，其他的数据为训练集。 Timestamp 可以支持逻辑比大小，所以我们就可以通过以下代码划分数据集。

split_dt = pd.Timestamp('2017-01-01 00:00:00')
train_df = raw_df[raw_df['ds']< split_dt]
test_df = raw_df[raw_df['ds']>= split_dt]

prophet建模

prophet 遵从sklearn 模型的API接口规范，也就是建模经典三部曲：创建模型--> fit训练数据--> predict 测试集。

代码还真的只有三行！！

model = Prophet()
model= model.fit(train_df)
y_hat  = model.predict(test_df)

本来想调用matplotlib 画图试试，结果看到prophet已经自带plot函数。

model.plot(y_hat)

结果分析

从人眼可见的角度来判断，预测的有点模样，有周期趋势，均值似乎也差不多。我们来看一下，误差多大，这里采用mean_absolute_error来衡量，并且计算误差占均值的比例。

print(mean_absolute_error(test_df['y'].values,y_hat['yhat']))
print(mean_absolute_error(test_df['y'].values,y_hat['yhat'])/test_df['y'].values.mean())


1782.2521239509801
0.12056939407796417

平均误差为1782MW,误差约等于观测对象水平的12%。 那这个误差是大还是小呢？模型是好还不是不好呢？ 我觉得单看准确度，肯定不是好的模型。但是问题出在哪里呢？是prophet 模型的上限就是这样的准确度？还是我们没有好好调整prophet的参数？

没有对比就没有伤害！！下一篇，我们可以对比一下我们的快速prophet 模型与网上的大牛的结果对比。

bonus：一点知识

到此我们就已经预测了最后一年的时间趋势，但是我们对prophet还是不甚了解。 prophet 核心就是时间序列的分解技术，它可以从原始序列中剥离每日周期序列，每周周期序列，每年周期序列，以及大趋势。

我们可以看看plot这次分解的结果。似乎大的趋势是一条水平线，也就是常值。每周的序列周期不是特别明显，但是每年的周期比较明显，周期约为半年，幅值高达2000左右。

fig2 = model.plot_components(y_hat)

总结

本文中我们采用prophet模型，对能源消耗的数据集的进行了预测，其中预测一年的平均误差为1782 左右。 本文还涵盖了使用prophet建模的一些辅助事项：

  • 安装过程如何避免坑
  • 数据如何准备可以被prophet 作为输入
  • plot prophet 分解的趋势和周期 由于本文讲求的快速建模，所以没有涵盖prophet的理论以及背景知识，网上这方面的资料很多，首推官方手册。
avatar
文章
阅读
粉丝