Python时间序列数据分析--以示例说明

时间  2021-02-16
标签 html python git 程序员 github 算法 微信 机器学习 函数 post 栏目 Python 繁體版
原文   原文链接


向AI转型的程序员都关注了这个号👇👇👇html

人工智能大数据与深度学习  公众号： datayxpython


本文的内容主要来源于博客：
https://www.analyticsvidhya.com/blog/2016/02/time-series-forecasting-codes-python/ 英文不错的读者能够前去阅读原文。git



在阅读本文以前 ，推荐先阅读：


时间序列预测之--ARIMA模型

http://www.cnblogs.com/bradleon/p/6827109.html


本文主要分为四个部分：程序员

  1. 用pandas处理时序数据github

  2. 怎样检查时序数据的稳定性算法

  3. 怎样让时序数据具备稳定性微信

  4. 时序数据的预测机器学习

1. 用pandas导入和处理时序数据

第一步：导入经常使用的库函数

import pandas as pd
import numpy as np
import matplotlib.pylab as plt
from matplotlib.pylab 
import rcParams
#rcParams设定好画布的大小

rcParams['figure.figsize'] = 15, 6

第二步：导入时序数据
数据文件可在github：post

关注微信公众号 datayx  而后 回复 “时序” 便可获取。


data = pd.read_csv(path+"AirPassengers.csv")
print data.head()
print '\n Data types:'
print data.dtypes


运行结果以下：数据包括每月对应的passenger的数目。
能够看到data已是一个DataFrame，包含两列Month和#Passengers，其中Month的类型是object，而index是0,1,2...


第三步：处理时序数据
咱们须要将Month的类型变为datetime，同时做为index。

dateparse = lambda dates: pd.datetime.strptime(dates, '%Y-%m')
#---其中parse_dates 代表选择数据中的哪一个column做为date-time信息，

#---index_col 告诉pandas以哪一个column做为 index

#--- date_parser 使用一个function(本文用lambda表达式代替)，使一个string转换为一个datetime变量

data = pd.read_csv('AirPassengers.csv', parse_dates=['Month'], index_col='Month',date_parser=dateparse)
print data.head()
print data.index


结果以下：能够看到data的index已经变成datetime类型的Month了。

2.怎样检查时序数据的稳定性(Stationarity)

由于ARIMA模型要求数据是稳定的，因此这一步相当重要。

1. 判断数据是稳定的常基于对于时间是常量的几个统计量：

  1. 常量的均值

  2. 常量的方差

  3. 与时间独立的自协方差

用图像说明以下：

  1. 均值

    X是时序数据的值，t是时间。能够看到左图，数据的均值对于时间轴来讲是常量，即数据的均值不是时间的函数,全部它是稳定的；右图随着时间的推移，数据的值总体趋势是增长的，全部均值是时间的函数，数据具备趋势，因此是非稳定的。

  2. 方差

    能够看到左图，数据的方差对于时间是常量，即数据的值域围绕着均值上下波动的振幅是固定的，因此左图数据是稳定的。而右图，数据的振幅在不一样时间点不一样，因此方差对于时间不是独立的，数据是非稳定的。可是左、右图的均值是一致的。


  3. 自协方差

    一个时序数据的自协方差，就是它在不一样两个时刻i,j的值的协方差。能够看到左图的自协方差于时间无关；而右图，随着时间的不一样，数据的波动频率明显不一样，致使它i，j取值不一样，就会获得不一样的协方差，所以是非稳定的。虽然右图在均值和方差上都是与时间无关的，但还是非稳定数据。


2. python判断时序数据稳定性

有两种方法：
1.Rolling statistic-- 即每一个时间段内的平均的数据均值和标准差状况。

  1. Dickey-Fuller Test -- 这个比较复杂，大体意思就是在必定置信水平下，对于时序数据假设 Null hypothesis: 非稳定。
    if 经过检验值(statistic)< 临界值(critical value)，则拒绝null hypothesis，即数据是稳定的；反之则是非稳定的。

from statsmodels.tsa.stattools import adfuller
def test_stationarity(timeseries):
    
    #这里以一年为一个窗口，每个时间t的值由它前面12个月（包括本身）的均值代替，标准差同理。
    rolmean = pd.rolling_mean(timeseries,window=12)
    rolstd = pd.rolling_std(timeseries, window=12)    
    #plot rolling statistics:
    fig = plt.figure()
    fig.add_subplot()
    orig = plt.plot(timeseries, color = 'blue',label='Original')
    mean = plt.plot(rolmean , color = 'red',label = 'rolling mean')
    std = plt.plot(rolstd, color = 'black', label= 'Rolling standard deviation')
    
    plt.legend(loc = 'best')
    plt.title('Rolling Mean & Standard Deviation')
    plt.show(block=False)    
    
    #Dickey-Fuller test:
    
    print 'Results of Dickey-Fuller Test:'
    dftest = adfuller(timeseries,autolag = 'AIC')    #dftest的输出前一项依次为检测值，p值，滞后数，使用的观测数，各个置信度下的临界值
    dfoutput = pd.Series(dftest[0:4],index = ['Test Statistic','p-value','#Lags Used','Number of Observations Used'])    
for key,value in dftest[4].items():
        dfoutput['Critical value (%s)' %key] = value    
    print dfoutput
    
ts = data['#Passengers']
test_stationarity(ts)


结果以下：


能够看到，数据的rolling均值/标准差具备愈来愈大的趋势，是不稳定的。
且DF-test能够明确的指出，在任何置信度下，数据都不是稳定的。


3. 让时序数据变成稳定的方法

让数据变得不稳定的缘由主要有俩：

  1. 趋势（trend）-数据随着时间变化。好比说升高或者下降。

  2. 季节性(seasonality)-数据在特定的时间段内变更。好比说节假日，或者活动致使数据的异常。

因为原数据值域范围比较大，为了缩小值域，同时保留其余信息，经常使用的方法是对数化，取log。

ts_log = np.log(ts)

检测和去除趋势
一般有三种方法：

    • 聚合 : 将时间轴缩短，以一段时间内星期/月/年的均值做为数据值。使不一样时间段内的值差距缩小。

    • 平滑： 以一个滑动窗口内的均值代替原来的值，为了使值之间的差距缩小

    • 多项式过滤：用一个回归模型来拟合现有数据，使得数据更平滑。

本文主要使用平滑方法

Moving Average--移动平均

moving_avg = pd.rolling_mean(ts_log,12)
plt.plot(ts_log ,color = 'blue')
plt.plot(moving_avg, color='red')



能够看出moving_average要比原值平滑许多。

而后做差：

ts_log_moving_avg_diff = ts_log-moving_avg
ts_log_moving_avg_diff.dropna(inplace = True)
test_stationarity(ts_log_moving_avg_diff)



能够看到，作了处理以后的数据基本上没有了随时间变化的趋势，DFtest的结果告诉咱们在95%的置信度下，数据是稳定的。

上面的方法是将全部的时间平等看待，而在许多状况下，能够认为越近的时刻越重要。因此引入指数加权移动平均-- Exponentially-weighted moving average.（pandas中经过ewma()函数提供了此功能。）

# halflife的值决定了衰减因子alpha：  alpha = 1 - exp(log(0.5) / halflife)expweighted_avg = pd.ewma(ts_log,halflife=12)
ts_log_ewma_diff = ts_log - expweighted_avg
test_stationarity(ts_log_ewma_diff)



能够看到相比普通的Moving Average，新的数据平均标准差更小了。并且DFtest能够获得结论：数据在99%的置信度上是稳定的。

  1. 检测和去除季节性
    有两种方法：

    • 1 差分化： 以特定滞后数目的时刻的值的做差

    • 2 分解： 对趋势和季节性分别建模在移除它们

Differencing--差分

ts_log_diff = ts_log - ts_log.shift()
ts_log_diff.dropna(inplace=True)
test_stationarity(ts_log_diff)


如图，能够看出相比MA方法，Differencing方法处理后的数据的均值和方差的在时间轴上的振幅明显缩小了。DFtest的结论是在90%的置信度下，数据是稳定的。


3.Decomposing-分解

#分解(decomposing) 能够用来把时序数据中的趋势和周期性数据都分离出来:

from statsmodels.tsa.seasonal import seasonal_decompose
def decompose(timeseries):
    
    # 返回包含三个部分 trend（趋势部分） ， seasonal（季节性部分） 和residual (残留部分)
    decomposition = seasonal_decompose(timeseries)
    
    trend = decomposition.trend
    seasonal = decomposition.seasonal
    residual = decomposition.resid
    
    plt.subplot(411)
    plt.plot(ts_log, label='Original')
    plt.legend(loc='best')
    plt.subplot(412)
    plt.plot(trend, label='Trend')
    plt.legend(loc='best')
    plt.subplot(413)
    plt.plot(seasonal,label='Seasonality')
    plt.legend(loc='best')
    plt.subplot(414)
    plt.plot(residual, label='Residuals')
    plt.legend(loc='best')
    plt.tight_layout()    
    return trend , seasonal, residual


如图能够明显的看到，将original数据 拆分红了三份。Trend数据具备明显的趋势性，Seasonality数据具备明显的周期性，Residuals是剩余的部分，能够认为是去除了趋势和季节性数据以后，稳定的数据，是咱们所须要的。

#消除了trend 和seasonal以后，只对residual部分做为想要的时序数据进行处理trend , seasonal, residual = decompose(ts_log)
residual.dropna(inplace=True)
test_stationarity(residual)



如图所示，数据的均值和方差趋于常数，几乎无波动(看上去比以前的陡峭，可是要注意他的值域只有[-0.05,0.05]之间)，因此直观上能够认为是稳定的数据。另外DFtest的结果显示，Statistic值原小于1%时的Critical value，因此在99%的置信度下，数据是稳定的。


4. 对时序数据进行预测

假设通过处理，已经获得了稳定时序数据。接下来，咱们使用ARIMA模型
对数据已经预测。ARIMA的介绍能够见本目录下的另外一篇文章。

step1： 经过ACF,PACF进行ARIMA（p，d，q）的p，q参数估计

由前文Differencing部分已知，一阶差分后数据已经稳定，因此d=1。
因此用一阶差分化的ts_log_diff = ts_log - ts_log.shift() 做为输入。
等价于