时间序列预测:探索性数据分析实用指南

2024年05月13日 由 alex 发表 82 0

简介

时间序列分析无疑是数据科学和机器学习领域最广泛的主题之一:无论是预测金融事件、能源消耗、产品销售还是股市趋势,这一领域一直是企业非常感兴趣的领域。


显然,数据可用性的大幅提高,加上机器学习模型的不断进步,使得这一话题如今更加引人关注。除了源自统计的传统预测方法(如回归模型、ARIMA 模型、指数平滑法)外,与机器学习(如基于树的模型)和深度学习(如 LSTM 网络、CNN、基于变换器的模型)相关的技术也已出现了一段时间。


尽管这些技术之间存在巨大差异,但无论采用何种模型,都必须完成一个初步步骤: 探索性数据分析。


在统计学中,探索性数据分析(EDA)是一门对数据进行分析和可视化的学科,目的是总结数据的主要特征并从中获取相关信息。这在数据科学领域相当重要,因为它可以为另一个重要步骤:特征工程奠定基础。也就是从数据集中创建、转换和提取特征,从而使模型发挥最大作用的实践。


因此,本文的目的是定义一个清晰的探索性数据分析模板,该模板侧重于时间序列,可以总结和突出数据集的最重要特征。为此,我们将使用一些常用的 Python 库,如 Pandas、Seaborn 和 Statsmodel。


数据

首先让我们定义一下数据集:在本文中,我们将使用 Kaggle 的每小时能源消耗数据。该数据集与 PJM 每小时能源消耗数据有关,PJM 是美国的一个区域输电组织,为特拉华州、伊利诺伊州、印第安纳州、肯塔基州、马里兰州、密歇根州、新泽西州、北卡罗来纳州、俄亥俄州、宾夕法尼亚州、田纳西州、弗吉尼亚州、西弗吉尼亚州和哥伦比亚特区提供电力服务。


每小时电力消耗数据来自 PJM 网站,单位为兆瓦 (MW)。


探索性数据分析

现在让我们来定义一下在处理时间序列时需要进行的最重要的分析。


当然,最重要的一点是绘制数据图:图表可以突出许多特征,如模式、异常观测、随时间的变化以及变量之间的关系。如前所述,从这些图表中得出的见解必须尽可能纳入预测模型。此外,一些数学工具,如描述性统计和时间序列分解,也会非常有用。


我在本文中提出的 EDA 包括六个步骤: 描述性统计、时间图、季节图、箱形图、时间序列分解、滞后分析。


1. 描述性统计

描述性统计是一种总结性统计,用于定量描述或总结结构化数据集合的特征。


常用于描述数据集的指标有:中心倾向度量(如平均值、中位数)、离散度量(如范围、标准差)和位置度量(如百分位数、四分位数)。所有这些都可以用所谓的五数汇总来概括,其中包括:分布的最小值、第一四分位数 (Q1)、中位数或第二四分位数 (Q2)、第三四分位数 (Q3) 和最大值。


在 Python 中,可以使用 Pandas 中众所周知的 describe 方法轻松获取这些信息:


import pandas as pd
# Loading and preprocessing steps
df = pd.read_csv('../input/hourly-energy-consumption/PJME_hourly.csv')
df = df.set_index('Datetime')
df.index = pd.to_datetime(df.index)
df.describe()


2


2. 时间图

首先要绘制的明显图形是时间图。也就是说,观察结果与观察时间相对应,连续的观察结果用线连接起来。


在 Python 中,我们可以使用 Pandas 和 Matplotlib:


import matplotlib.pyplot as plt
 
# Set pyplot style
plt.style.use("seaborn")
# Plot
df['PJME_MW'].plot(title='PJME - Time Plot', figsize=(10,6))
plt.ylabel('Consumption [MW]')
plt.xlabel('Date')


3


这幅图已经提供了一些信息:


  1. 正如我们所预料的那样,这种模式显示出每年的季节性。
  2. 如果只看一年,似乎会发现更多的规律。由于耗电量较大,冬季和夏季的耗电量很可能各有一个高峰。
  3. 多年来,该系列并没有表现出明显的增减趋势,平均消耗量保持稳定。
  4. 在 2023 年前后有一个异常值,可能应在实施模型时加以估算。


3. 季节图

季节图从根本上说是一种时间图,它将数据与所属系列的各个 "季节 "相对应。


就能源消耗而言,我们通常可以获得每小时的数据,因此可能存在多种季节性:每年、每周、每天。在深入研究这些图表之前,让我们先在 Pandas 数据框中设置一些变量:


# Defining required fields
df['year'] = [x for x in df.index.year]
df['month'] = [x for x in df.index.month]
df = df.reset_index()
df['week'] = df['Datetime'].apply(lambda x:x.week)
df = df.set_index('Datetime')
df['hour'] = [x for x in df.index.hour]
df['day'] = [x for x in df.index.day_of_week]
df['day_str'] = [x.strftime('%a') for x in df.index]
df['year_month'] = [str(x.year) + '_' + str(x.month) for x in df.index]


3.1 季节性曲线图--年度消耗量

一个非常有趣的图示是按年和月分组的能源消耗图,它突出了每年的季节性,并能告诉我们这些年来的上升/下降趋势。


以下是 Python 代码:


import numpy as np
# Defining colors palette
np.random.seed(42)
df_plot = df[['month', 'year', 'PJME_MW']].dropna().groupby(['month', 'year']).mean()[['PJME_MW']].reset_index()
years = df_plot['year'].unique()
colors = np.random.choice(list(mpl.colors.XKCD_COLORS.keys()), len(years), replace=False)
# Plot
plt.figure(figsize=(16,12))
for i, y in enumerate(years):
    if i > 0:        
        plt.plot('month', 'PJME_MW', data=df_plot[df_plot['year'] == y], color=colors[i], label=y)
        if y == 2018:
            plt.text(df_plot.loc[df_plot.year==y, :].shape[0]+0.3, df_plot.loc[df_plot.year==y, 'PJME_MW'][-1:].values[0], y, fontsize=12, color=colors[i])
        else:
            plt.text(df_plot.loc[df_plot.year==y, :].shape[0]+0.1, df_plot.loc[df_plot.year==y, 'PJME_MW'][-1:].values[0], y, fontsize=12, color=colors[i])
# Setting labels
plt.gca().set(ylabel= 'PJME_MW', xlabel = 'Month')
plt.yticks(fontsize=12, alpha=.7)
plt.title("Seasonal Plot - Monthly Consumption", fontsize=20)
plt.ylabel('Consumption [MW]')
plt.xlabel('Month')
plt.show()


4


从图中可以看出,每年的消耗量实际上都有一个非常固定的模式:冬季消耗量明显增加,夏季消耗量达到峰值(由于供暖/制冷系统),而春季和秋季通常不需要供暖或制冷,因此消耗量为最小值。


此外,这幅图还告诉我们,不同年份的总体消耗量并没有明显的增减模式。


3.2 季节图--每周消耗量

另一种有用的曲线图是周曲线图,它描绘的是一周内不同月份的消耗量,也可以说明一年内每周消耗量是否在变化以及如何变化。


让我们看看如何用 Python 来计算:


# Defining colors palette
np.random.seed(42)
df_plot = df[['month', 'day_str', 'PJME_MW', 'day']].dropna().groupby(['day_str', 'month', 'day']).mean()[['PJME_MW']].reset_index()
df_plot = df_plot.sort_values(by='day', ascending=True)
months = df_plot['month'].unique()
colors = np.random.choice(list(mpl.colors.XKCD_COLORS.keys()), len(months), replace=False)
# Plot
plt.figure(figsize=(16,12))
for i, y in enumerate(months):
    if i > 0:        
        plt.plot('day_str', 'PJME_MW', data=df_plot[df_plot['month'] == y], color=colors[i], label=y)
        if y == 2018:
            plt.text(df_plot.loc[df_plot.month==y, :].shape[0]-.9, df_plot.loc[df_plot.month==y, 'PJME_MW'][-1:].values[0], y, fontsize=12, color=colors[i])
        else:
            plt.text(df_plot.loc[df_plot.month==y, :].shape[0]-.9, df_plot.loc[df_plot.month==y, 'PJME_MW'][-1:].values[0], y, fontsize=12, color=colors[i])

# Setting Labels
plt.gca().set(ylabel= 'PJME_MW', xlabel = 'Month')
plt.yticks(fontsize=12, alpha=.7)
plt.title("Seasonal Plot - Weekly Consumption", fontsize=20)
plt.ylabel('Consumption [MW]')
plt.xlabel('Month')
plt.show()


5


3.3 季节性曲线图--日消耗量

最后,我要展示的最后一个季节图是日消耗量图。正如你所猜测的那样,它代表了一天中消费量的变化情况。在这种情况下,数据首先按星期分组,然后取平均值汇总。


下面是代码


import seaborn as sns
# Defining the dataframe
df_plot = df[['hour', 'day_str', 'PJME_MW']].dropna().groupby(['hour', 'day_str']).mean()[['PJME_MW']].reset_index()
# Plot using Seaborn
plt.figure(figsize=(10,8))
sns.lineplot(data = df_plot, x='hour', y='PJME_MW', hue='day_str', legend=True)
plt.locator_params(axis='x', nbins=24)
plt.title("Seasonal Plot - Daily Consumption", fontsize=20)
plt.ylabel('Consumption [MW]')
plt.xlabel('Hour')
plt.legend()


6


通常情况下,这种曲线图会显示出一种非常典型的模式,有人称之为 "M 型曲线",因为白天的消耗量似乎呈 "M "型。有时这种模式很明显,有时则不明显(比如本例)。


不过,这种曲线图通常会在一天的中间(上午 10 点到下午 2 点)出现一个相对的峰值,然后是一个相对的峰值(下午 2 点到下午 6 点)和另一个峰值(下午 6 点到晚上 8 点)。最后,它还显示了周末和其他日子的消耗量差异。


3.4 季节图--特征工程

现在我们来看看如何将这些信息用于特征工程。假设我们正在使用一些需要高质量特征的 ML 模型(如 ARIMA 模型或基于树的模型)。


这些是来自季节图的主要证据:


  1. 年度消费在不同年份变化不大:这表明可以使用来自滞后变量或外生变量的年度季节性特征。
  2. 各月的周消费遵循相同的模式:这表明可以使用来自滞后变量或外生变量的周特征。
  3. 日常消费与平日和周末不同:这表明可以使用分类特征来识别平日和非平日。


4. 方框图

方框图是确定数据分布情况的有效方法。简而言之,方框图描述了百分位数(代表分布的第 1 个四分位数(Q1)、第 2 个四分位数(Q2/中位数)和第 3 个四分位数(Q3))和边线(代表数据的范围)。超出晶须的每一个值都可以被视为离群值,更深入地说,晶须的计算方法通常是:


7


4.1 箱形图 - 总消耗量

我们首先来计算总消耗量的箱形图,这可以通过 Seaborn 轻松完成:


plt.figure(figsize=(8,5))8,5))
sns.boxplot(data=df, x='PJME_MW')
plt.xlabel('Consumption [MW]')
plt.title(f'Boxplot - Consumption Distribution');


8


尽管这幅图看起来信息量不大,但它告诉我们,我们正在处理的是一个类似高斯分布的 分布,其尾部更偏向右侧。


4.2 箱形图--日月分布

日/月箱形图是一个非常有趣的图。通过创建一个 "日-月 "变量,并根据该变量对消耗量进行分组,可以得到该图。以下是代码,仅涉及 2017 年:


df['year'] = [x for x in df.index.year]'year'] = [x for x in df.index.year]
df['month'] = [x for x in df.index.month]
df['year_month'] = [str(x.year) + '_' + str(x.month) for x in df.index]
df_plot = df[df['year'] >= 2017].reset_index().sort_values(by='Datetime').set_index('Datetime')
plt.title(f'Boxplot Year Month Distribution');
plt.xticks(rotation=90)
plt.xlabel('Year Month')
plt.ylabel('MW')
sns.boxplot(x='year_month', y='PJME_MW', data=df_plot)
plt.ylabel('Consumption [MW]')
plt.xlabel('Year Month')


9


可以看出,夏/冬两季(即气温达到峰值时)的消费量不确定性较小,而春/秋两季(即气温变化较大时)的消费量较为分散。最后,2018 年夏季的消费量高于 2017 年,这可能是由于夏季较为温暖的缘故。在进行特征工程设计时,请记住将温度曲线(如果有的话)包括在内,因为它可能可以用作外生变量。


4.3 方框图--日分布

另一种有用的曲线图是指一周内的消耗量分布,这与周消耗量季节曲线图类似。


df_plot = df[['day_str', 'day', 'PJME_MW']].sort_values(by='day')'day_str', 'day', 'PJME_MW']].sort_values(by='day')
plt.title(f'Boxplot Day Distribution')
plt.xlabel('Day of week')
plt.ylabel('MW')
sns.boxplot(x='day_str', y='PJME_MW', data=df_plot)
plt.ylabel('Consumption [MW]')
plt.xlabel('Day of week')


10


如前所述,周末的消费量明显较低。总之,有几个异常值表明,"星期 "等日历特征肯定是有用的,但不能完全解释这一系列现象。


4.4 箱形图--小时分布

最后让我们来看看小时分布箱形图。它与每日消费季节图类似,因为它提供了消费在一天中的分布情况。代码如下


plt.title(f'Boxplot Hour Distribution');f'Boxplot Hour Distribution');
plt.xlabel('Hour')
plt.ylabel('MW')
sns.boxplot(x='hour', y='PJME_MW', data=df)
plt.ylabel('Consumption [MW]')
plt.xlabel('Hour')


11


请注意,之前看到的 "M "形现在变得更粗了。此外,还有很多离群值,这说明数据不仅依赖于每日的季节性(例如,今天上午 12 点的消耗量与昨天上午 12 点的消耗量相似),还依赖于其他因素,可能是温度或湿度等外生气候特征。


5. 时间序列分解

如前所述,时间序列数据可以表现出多种模式。通常情况下,将一个时间序列拆分成几个部分是有帮助的,每个部分代表一个基本模式类别。


我们可以把时间序列看作由三个部分组成:趋势部分、季节部分和剩余部分(包含时间序列中的任何其他部分)。对于某些时间序列(如能源消耗序列),可能存在不止一个季节性成分,分别对应于不同的季节性时期(日、周、月、年)。


分解方法主要有两种:加法分解和乘法分解。


对于加法分解,我们将一个序列(?)表示为季节成分(?)、趋势(?)和余数(?)的总和:


12


同样,乘法分解可以写成


13


一般来说,加法分解最能代表方差恒定的序列,而乘法分解最适合方差非平稳的时间序列。


在 Python 中,使用 Statsmodel 库可以轻松实现时间序列分解:


df_plot = df[df['year'] == 2017].reset_index()'year'] == 2017].reset_index()
df_plot = df_plot.drop_duplicates(subset=['Datetime']).sort_values(by='Datetime')
df_plot = df_plot.set_index('Datetime')
df_plot['PJME_MW - Multiplicative Decompose'] = df_plot['PJME_MW']
df_plot['PJME_MW - Additive Decompose'] = df_plot['PJME_MW']
# Additive Decomposition
result_add = seasonal_decompose(df_plot['PJME_MW - Additive Decompose'], model='additive', period=24*7)
# Multiplicative Decomposition 
result_mul = seasonal_decompose(df_plot['PJME_MW - Multiplicative Decompose'], model='multiplicative', period=24*7)
# Plot
result_add.plot().suptitle('', fontsize=22)
plt.xticks(rotation=45)
result_mul.plot().suptitle('', fontsize=22)
plt.xticks(rotation=45)
plt.show()


14


15


上图是 2017 年的数据。在这两种情况下,我们都可以看到趋势有几个局部峰值,夏季数值较高。从季节性成分中,我们可以看到该序列实际上有几个周期性,这幅图更多地突出了周周期性,但如果我们关注同年的某个特定月份(1 月),日周期性也会出现:


df_plot = df[(df['year'] == 2017)].reset_index()
df_plot = df_plot[df_plot['month'] == 1]
df_plot['PJME_MW - Multiplicative Decompose'] = df_plot['PJME_MW']
df_plot['PJME_MW - Additive Decompose'] = df_plot['PJME_MW']
df_plot = df_plot.drop_duplicates(subset=['Datetime']).sort_values(by='Datetime')
df_plot = df_plot.set_index('Datetime')
# Additive Decomposition
result_add = seasonal_decompose(df_plot['PJME_MW - Additive Decompose'], model='additive', period=24*7)
# Multiplicative Decomposition 
result_mul = seasonal_decompose(df_plot['PJME_MW - Multiplicative Decompose'], model='multiplicative', period=24*7)
# Plot
result_add.plot().suptitle('', fontsize=22)
plt.xticks(rotation=45)
result_mul.plot().suptitle('', fontsize=22)
plt.xticks(rotation=45)
plt.show()


16


17


6. 滞后分析

在时间序列预测中,滞后期就是序列的过去值。例如,对于日序列,第一个滞后值指的是该序列前一天的值,第二个滞后值指的是前一天的值,以此类推。


滞后分析的基础是计算序列与序列本身的滞后版本之间的相关性,这也称为自相关性。对于一个 k 个滞后版本的序列,我们将自相关系数定义为


18


其中,y 表示序列的平均值,k 表示滞后期。


自相关系数构成了序列的自相关函数(ACF),即自相关系数与滞后期数的关系图。


当数据具有趋势性时,小滞后期的自相关系数通常较大且为正,因为时间上接近的观测值在数值上也接近。当数据具有季节性时,与季节性滞后期(和季节性周期的倍数)相对应的自相关值会比其他滞后期大。同时具有趋势和季节性的数据将显示这些效应的组合。


实际上,更有用的函数是部分自相关函数(PACF)。它与 ACF 类似,只是只显示两个滞后期之间的直接自相关。例如,滞后 3 的部分自相关指的是滞后 1 和滞后 2 无法解释的唯一相关性。换句话说,部分相关指的是某个滞后期对当前时间值的直接影响。


在开始 Python 代码之前,需要强调的是,如果序列是静态的,自相关系数会更明显,因此最好先对序列进行分化,以稳定信号。


下面是绘制一天中不同时段 PACF 的代码:


from statsmodels.graphics.tsaplots import plot_pacf
actual = df['PJME_MW']
hours = range(0, 24, 4)
for hour in hours:
    plot_pacf(actual[actual.index.hour == hour].diff().dropna(), lags=30, alpha=0.01)
    plt.title(f'PACF - h = {hour}')
    plt.ylabel('Correlation')
    plt.xlabel('Lags')
    plt.show()


19


20


21


22


23


24


正如你所看到的,PACF 简单来说就是绘制不同滞后期的皮尔逊部分自相关系数。当然,非滞后序列本身显示出完美的自相关性,因此滞后 0 总是 1。蓝色区域代表置信区间:如果滞后期超过该区域,则在统计上具有显著性,我们可以断言其具有重要意义。


6.1 滞后分析--特征工程

滞后分析是对时间序列特征工程影响最大的研究之一。如前所述,相关性高的滞后期是序列的重要滞后期,因此应加以考虑。


一种广泛使用的特征工程技术是按小时划分数据集。也就是说,将数据分成 24 个子集,每个子集指一天中的一个小时。这样做的效果是使信号正则化和平滑化,从而使预测更加简单。


然后对每个子集进行特征设计、训练和微调。最终的预测将综合这 24 个模型的结果。也就是说,每个小时模型都有其特殊性,其中大部分都会考虑到重要的滞后因素。


在继续讨论之前,我们先来定义一下滞后分析中的两种滞后情况:


  1. 自回归滞后期:接近滞后期 0 的滞后期,我们预期其值较高(最近的滞后期更有可能预测现值)。自回归滞后期代表了序列的趋势程度。
  2. 季节滞后期:滞后期指的是季节性时期。当按小时分割数据时,它们通常代表每周的季节性。


现在我们来讨论上面打印的 PACF 图。


夜间时间

夜间时间(0,4)的消费更依赖于自回归滞后期而不是周滞后期,因为最重要的都集中在前五个滞后期。7、14、21、28 等季节性时段似乎不太重要,这建议我们在进行特征工程时特别关注滞后期 1 至 5。


日间时段

日间时段(8、12、16、20)的消耗量既表现出自回归滞后,也表现出季节性滞后。这一点在 8 和 12 小时尤为明显,因为这两个小时的消费量特别高,而在临近夜间时,季节性滞后就变得不那么重要了。对于这些子集,我们还应该包括季节滞后和自回归滞后。


结论

本文旨在为时间序列预测提供一个全面的探索性数据分析模板。


在任何类型的数据科学研究中,EDA 都是一个基本步骤,因为它可以了解数据的性质和特殊性,并为特征工程奠定基础,而特征工程反过来又可以显著提高模型性能。



文章来源:https://medium.com/towards-data-science/time-series-forecasting-a-practical-guide-to-exploratory-data-analysis-a101dc5f85b1
欢迎关注ATYUN官方公众号
商务合作及内容投稿请联系邮箱:bd@atyun.com
评论 登录
热门职位
Maluuba
20000~40000/月
Cisco
25000~30000/月 深圳市
PilotAILabs
30000~60000/年 深圳市
写评论取消
回复取消