热门标签 | HotTags
当前位置:  开发笔记 > 编程语言 > 正文

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

本文的内容主要来源于博客:本人做了适当的注释和补充。https:www.analyticsvidhya.comblog201602time-series-forec


本文的内容主要来源于博客:本人做了适当的注释和补充。


https://www.analyticsvidhya.com/blog/2016/02/time-series-forecasting-codes-python/ 英文不错的读者可以前去阅读原文。

在阅读本文之前 ,推荐先阅读:http://www.cnblogs.com/bradleon/p/6827109.html


导读

本文主要分为四个部分:


  1. 用pandas处理时序数据
  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:
http://github.com/aarshayj/Analytics_Vidhya/tree/master/Articles/Time_Series_Analysis 中下载


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...
filelist

第三步:处理时序数据
我们需要将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了。
filelist
filelist


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

因为ARIMA模型要求数据是稳定的,所以这一步至关重要。


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


  1. 常量的均值
  2. 常量的方差
  3. 与时间独立的自协方差

用图像说明如下:


  1. 均值
    filelist
    X是时序数据的值,t是时间。可以看到左图,数据的均值对于时间轴来说是常量,即数据的均值不是时间的函数,所有它是稳定的;右图随着时间的推移,数据的值整体趋势是增加的,所有均值是时间的函数,数据具有趋势,所以是非稳定的。
  2. 方差
    filelist
    可以看到左图,数据的方差对于时间是常量,即数据的值域围绕着均值上下波动的振幅是固定的,所以左图数据是稳定的。而右图,数据的振幅在不同时间点不同,所以方差对于时间不是独立的,数据是非稳定的。但是左、右图的均值是一致的。
  3. 自协方差
    filelist
    一个时序数据的自协方差,就是它在不同两个时刻i,j的值的协方差。可以看到左图的自协方差于时间无关;而右图,随着时间的不同,数据的波动频率明显不同,导致它i,j取值不同,就会得到不同的协方差,因此是非稳定的。虽然右图在均值和方差上都是与时间无关的,但仍是非稳定数据。

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

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


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

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

结果如下&#xff1a;


可以看到&#xff0c;数据的rolling均值/标准差具有越来越大的趋势&#xff0c;是不稳定的。
且DF-test可以明确的指出&#xff0c;在任何置信度下&#xff0c;数据都不是稳定的。


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

让数据变得不稳定的原因主要有俩&#xff1a;


  1. 趋势&#xff08;trend&#xff09;-数据随着时间变化。比如说升高或者降低。
  2. 季节性(seasonality)-数据在特定的时间段内变动。比如说节假日&#xff0c;或者活动导致数据的异常。

由于原数据值域范围比较大&#xff0c;为了缩小值域&#xff0c;同时保留其他信息&#xff0c;常用的方法是对数化&#xff0c;取log。


ts_log &#61; np.log(ts)


  1. 检测和去除趋势
    通常有三种方法&#xff1a;

    • 聚合 : 将时间轴缩短&#xff0c;以一段时间内星期/月/年的均值作为数据值。使不同时间段内的值差距缩小。
    • 平滑&#xff1a; 以一个滑动窗口内的均值代替原来的值&#xff0c;为了使值之间的差距缩小
    • 多项式过滤&#xff1a;用一个回归模型来拟合现有数据&#xff0c;使得数据更平滑。

本文主要使用平滑方法


Moving Average--移动平均


moving_avg &#61; pd.rolling_mean(ts_log,12)
plt.plot(ts_log ,color &#61; &#39;blue&#39;)
plt.plot(moving_avg, color&#61;&#39;red&#39;)


可以看出moving_average要比原值平滑许多。

然后作差&#xff1a;


ts_log_moving_avg_diff &#61; ts_log-moving_avg
ts_log_moving_avg_diff.dropna(inplace &#61; True)
test_stationarity(ts_log_moving_avg_diff)



可以看到&#xff0c;做了处理之后的数据基本上没有了随时间变化的趋势&#xff0c;DFtest的结果告诉我们在95%的置信度下&#xff0c;数据是稳定的。

上面的方法是将所有的时间平等看待&#xff0c;而在许多情况下&#xff0c;可以认为越近的时刻越重要。所以引入指数加权移动平均-- Exponentially-weighted moving average.&#xff08;pandas中通过ewma()函数提供了此功能。&#xff09;


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



可以看到相比普通的Moving Average&#xff0c;新的数据平均标准差更小了。而且DFtest可以得到结论&#xff1a;数据在99%的置信度上是稳定的。


  1. 检测和去除季节性
    有两种方法&#xff1a;

    • 1 差分化&#xff1a; 以特定滞后数目的时刻的值的作差
    • 2 分解&#xff1a; 对趋势和季节性分别建模在移除它们

Differencing--差分


ts_log_diff &#61; ts_log - ts_log.shift()
ts_log_diff.dropna(inplace&#61;True)
test_stationarity(ts_log_diff)


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

3.Decomposing-分解


#分解(decomposing) 可以用来把时序数据中的趋势和周期性数据都分离出来:
from statsmodels.tsa.seasonal import seasonal_decompose
def decompose(timeseries):# 返回包含三个部分 trend&#xff08;趋势部分&#xff09; &#xff0c; seasonal&#xff08;季节性部分&#xff09; 和residual (残留部分)decomposition &#61; seasonal_decompose(timeseries)trend &#61; decomposition.trendseasonal &#61; decomposition.seasonalresidual &#61; decomposition.residplt.subplot(411)plt.plot(ts_log, label&#61;&#39;Original&#39;)plt.legend(loc&#61;&#39;best&#39;)plt.subplot(412)plt.plot(trend, label&#61;&#39;Trend&#39;)plt.legend(loc&#61;&#39;best&#39;)plt.subplot(413)plt.plot(seasonal,label&#61;&#39;Seasonality&#39;)plt.legend(loc&#61;&#39;best&#39;)plt.subplot(414)plt.plot(residual, label&#61;&#39;Residuals&#39;)plt.legend(loc&#61;&#39;best&#39;)plt.tight_layout()return trend , seasonal, residual


如图可以明显的看到&#xff0c;将original数据 拆分成了三份。Trend数据具有明显的趋势性&#xff0c;Seasonality数据具有明显的周期性&#xff0c;Residuals是剩余的部分&#xff0c;可以认为是去除了趋势和季节性数据之后&#xff0c;稳定的数据&#xff0c;是我们所需要的。


#消除了trend 和seasonal之后&#xff0c;只对residual部分作为想要的时序数据进行处理
trend , seasonal, residual &#61; decompose(ts_log)
residual.dropna(inplace&#61;True)
test_stationarity(residual)



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


4. 对时序数据进行预测

假设经过处理&#xff0c;已经得到了稳定时序数据。接下来&#xff0c;我们使用ARIMA模型
对数据已经预测。ARIMA的介绍可以见本目录下的另一篇文章。

step1&#xff1a; 通过ACF,PACF进行ARIMA&#xff08;p&#xff0c;d&#xff0c;q&#xff09;的p&#xff0c;q参数估计

由前文Differencing部分已知&#xff0c;一阶差分后数据已经稳定&#xff0c;所以d&#61;1。
所以用一阶差分化的ts_log_diff &#61; ts_log - ts_log.shift() 作为输入。
等价于


yt&#61;YtYt1yt&#61;Yt−Yt−1
作为输入。

先画出ACF,PACF的图像,代码如下&#xff1a;


#ACF and PACF plots:
from statsmodels.tsa.stattools import acf, pacf
lag_acf &#61; acf(ts_log_diff, nlags&#61;20)
lag_pacf &#61; pacf(ts_log_diff, nlags&#61;20, method&#61;&#39;ols&#39;)
#Plot ACF:
plt.subplot(121)
plt.plot(lag_acf)
plt.axhline(y&#61;0,linestyle&#61;&#39;--&#39;,color&#61;&#39;gray&#39;)
plt.axhline(y&#61;-1.96/np.sqrt(len(ts_log_diff)),linestyle&#61;&#39;--&#39;,color&#61;&#39;gray&#39;)
plt.axhline(y&#61;1.96/np.sqrt(len(ts_log_diff)),linestyle&#61;&#39;--&#39;,color&#61;&#39;gray&#39;)
plt.title(&#39;Autocorrelation Function&#39;)#Plot PACF:
plt.subplot(122)
plt.plot(lag_pacf)
plt.axhline(y&#61;0,linestyle&#61;&#39;--&#39;,color&#61;&#39;gray&#39;)
plt.axhline(y&#61;-1.96/np.sqrt(len(ts_log_diff)),linestyle&#61;&#39;--&#39;,color&#61;&#39;gray&#39;)
plt.axhline(y&#61;1.96/np.sqrt(len(ts_log_diff)),linestyle&#61;&#39;--&#39;,color&#61;&#39;gray&#39;)
plt.title(&#39;Partial Autocorrelation Function&#39;)
plt.tight_layout()


图中&#xff0c;上下两条灰线之间是置信区间&#xff0c;p的值就是ACF第一次穿过上置信区间时的横轴值。q的值就是PACF第一次穿过上置信区间的横轴值。所以从图中可以得到p&#61;2&#xff0c;q&#61;2。

step2&#xff1a; 得到参数估计值p&#xff0c;d&#xff0c;q之后&#xff0c;生成模型ARIMA&#xff08;p&#xff0c;d&#xff0c;q&#xff09;
为了突出差别&#xff0c;用三种参数取值的三个模型作为对比。
模型1&#xff1a;AR模型(ARIMA(2,1,0))


from statsmodels.tsa.arima_model import ARIMA
model &#61; ARIMA(ts_log, order&#61;(2, 1, 0))
results_AR &#61; model.fit(disp&#61;-1)
plt.plot(ts_log_diff)
plt.plot(results_AR.fittedvalues, color&#61;&#39;red&#39;)
plt.title(&#39;RSS: %.4f&#39;% sum((results_AR.fittedvalues-ts_log_diff)**2))


图中&#xff0c;蓝线是输入值&#xff0c;红线是模型的拟合值&#xff0c;RSS的累计平方误差。

模型2&#xff1a;MA模型&#xff08;ARIMA&#xff08;0,1,2&#xff09;&#xff09;


model &#61; ARIMA(ts_log, order&#61;(0, 1, 2))
results_MA &#61; model.fit(disp&#61;-1)
plt.plot(ts_log_diff)
plt.plot(results_MA.fittedvalues, color&#61;&#39;red&#39;)
plt.title(&#39;RSS: %.4f&#39;% sum((results_MA.fittedvalues-ts_log_diff)**2))

模型3&#xff1a;ARIMA模型(ARIMA(2,1,2))


model &#61; ARIMA(ts_log, order&#61;(2, 1, 2))
results_ARIMA &#61; model.fit(disp&#61;-1)
plt.plot(ts_log_diff)
plt.plot(results_ARIMA.fittedvalues, color&#61;&#39;red&#39;)
plt.title(&#39;RSS: %.4f&#39;% sum((results_ARIMA.fittedvalues-ts_log_diff)**2))


由RSS&#xff0c;可知模型3--ARIMA&#xff08;2,1,2&#xff09;的拟合度最好&#xff0c;所以我们确定了最终的预测模型。

step3: 将模型代入原数据进行预测
因为上面的模型的拟合值是对原数据进行稳定化之后的输入数据的拟合&#xff0c;所以需要对拟合值进行相应处理的逆操作&#xff0c;使得它回到与原数据一致的尺度。



#ARIMA拟合的其实是一阶差分ts_log_diff&#xff0c;predictions_ARIMA_diff[i]是第i个月与i-1个月的ts_log的差值。
#由于差分化有一阶滞后&#xff0c;所以第一个月的数据是空的&#xff0c;
predictions_ARIMA_diff &#61; pd.Series(results_ARIMA.fittedvalues, copy&#61;True)
print predictions_ARIMA_diff.head()
#累加现有的diff&#xff0c;得到每个值与第一个月的差分&#xff08;同log底的情况下&#xff09;。
#即predictions_ARIMA_diff_cumsum[i] 是第i个月与第1个月的ts_log的差值。
predictions_ARIMA_diff_cumsum &#61; predictions_ARIMA_diff.cumsum()
#先ts_log_diff &#61;> ts_log&#61;>ts_log &#61;> ts
#先以ts_log的第一个值作为基数&#xff0c;复制给所有值&#xff0c;然后每个时刻的值累加与第一个月对应的差值(这样就解决了&#xff0c;第一个月diff数据为空的问题了)
#然后得到了predictions_ARIMA_log &#61;> predictions_ARIMA
predictions_ARIMA_log &#61; pd.Series(ts_log.ix[0], index&#61;ts_log.index)
predictions_ARIMA_log &#61; predictions_ARIMA_log.add(predictions_ARIMA_diff_cumsum,fill_value&#61;0)
predictions_ARIMA &#61; np.exp(predictions_ARIMA_log)
plt.figure()
plt.plot(ts)
plt.plot(predictions_ARIMA)
plt.title(&#39;RMSE: %.4f&#39;% np.sqrt(sum((predictions_ARIMA-ts)**2)/len(ts)))


5.总结

前面一篇文章&#xff0c;总结了ARIMA建模的步骤。
(1). 获取被观测系统时间序列数据&#xff1b;
(2). 对数据绘图&#xff0c;观测是否为平稳时间序列&#xff1b;对于非平稳时间序列要先进行d阶差分运算&#xff0c;化为平稳时间序列&#xff1b;
(3). 经过第二步处理&#xff0c;已经得到平稳时间序列。要对平稳时间序列分别求得其自相关系数ACF 和偏自相关系数PACF&#xff0c;通过对自相关图和偏自相关图的分析&#xff0c;得到最佳的阶层 p 和阶数 q

(4). 由以上得到的d、q、p&#xff0c;得到ARIMA模型。然后开始对得到的模型进行模型检验。



文章出处&#xff1a;http://www.cnblogs.com/bradleon/p/6832867.html


文章出处&#xff1a;http://www.cnblogs.com/bradleon/p/6832867.html

推荐阅读
  • 本文介绍了在Python3中如何使用选择文件对话框的格式打开和保存图片的方法。通过使用tkinter库中的filedialog模块的asksaveasfilename和askopenfilename函数,可以方便地选择要打开或保存的图片文件,并进行相关操作。具体的代码示例和操作步骤也被提供。 ... [详细]
  • 闭包一直是Java社区中争论不断的话题,很多语言都支持闭包这个语言特性,闭包定义了一个依赖于外部环境的自由变量的函数,这个函数能够访问外部环境的变量。本文以JavaScript的一个闭包为例,介绍了闭包的定义和特性。 ... [详细]
  • Allegro总结:1.防焊层(SolderMask):又称绿油层,PCB非布线层,用于制成丝网印板,将不需要焊接的地方涂上防焊剂.在防焊层上预留的焊盘大小要比实际的焊盘大一些,其差值一般 ... [详细]
  • 安卓select模态框样式改变_微软Office风格的多端(Web、安卓、iOS)组件库——Fabric UI...
    介绍FabricUI是微软开源的一套Office风格的多端组件库,共有三套针对性的组件,分别适用于web、android以及iOS,Fab ... [详细]
  • 【shell】网络处理:判断IP是否在网段、两个ip是否同网段、IP地址范围、网段包含关系
    本文介绍了使用shell脚本判断IP是否在同一网段、判断IP地址是否在某个范围内、计算IP地址范围、判断网段之间的包含关系的方法和原理。通过对IP和掩码进行与计算,可以判断两个IP是否在同一网段。同时,还提供了一段用于验证IP地址的正则表达式和判断特殊IP地址的方法。 ... [详细]
  • 本文介绍了如何在Mac上使用Pillow库加载不同于默认字体和大小的字体,并提供了一个简单的示例代码。通过该示例,读者可以了解如何在Python中使用Pillow库来写入不同字体的文本。同时,本文也解决了在Mac上使用Pillow库加载字体时可能遇到的问题。读者可以根据本文提供的示例代码,轻松实现在Mac上使用Pillow库加载不同字体的功能。 ... [详细]
  • 本文介绍了利用ARMA模型对平稳非白噪声序列进行建模的步骤及代码实现。首先对观察值序列进行样本自相关系数和样本偏自相关系数的计算,然后根据这些系数的性质选择适当的ARMA模型进行拟合,并估计模型中的位置参数。接着进行模型的有效性检验,如果不通过则重新选择模型再拟合,如果通过则进行模型优化。最后利用拟合模型预测序列的未来走势。文章还介绍了绘制时序图、平稳性检验、白噪声检验、确定ARMA阶数和预测未来走势的代码实现。 ... [详细]
  • 动量|收益率_基于MT策略的实战分析
    篇首语:本文由编程笔记#小编为大家整理,主要介绍了基于MT策略的实战分析相关的知识,希望对你有一定的参考价值。基于MT策略的实战分析 ... [详细]
  • 在Android开发中,使用Picasso库可以实现对网络图片的等比例缩放。本文介绍了使用Picasso库进行图片缩放的方法,并提供了具体的代码实现。通过获取图片的宽高,计算目标宽度和高度,并创建新图实现等比例缩放。 ... [详细]
  • 云原生边缘计算之KubeEdge简介及功能特点
    本文介绍了云原生边缘计算中的KubeEdge系统,该系统是一个开源系统,用于将容器化应用程序编排功能扩展到Edge的主机。它基于Kubernetes构建,并为网络应用程序提供基础架构支持。同时,KubeEdge具有离线模式、基于Kubernetes的节点、群集、应用程序和设备管理、资源优化等特点。此外,KubeEdge还支持跨平台工作,在私有、公共和混合云中都可以运行。同时,KubeEdge还提供数据管理和数据分析管道引擎的支持。最后,本文还介绍了KubeEdge系统生成证书的方法。 ... [详细]
  • 目录实现效果:实现环境实现方法一:基本思路主要代码JavaScript代码总结方法二主要代码总结方法三基本思路主要代码JavaScriptHTML总结实 ... [详细]
  • [译]技术公司十年经验的职场生涯回顾
    本文是一位在技术公司工作十年的职场人士对自己职业生涯的总结回顾。她的职业规划与众不同,令人深思又有趣。其中涉及到的内容有机器学习、创新创业以及引用了女性主义者在TED演讲中的部分讲义。文章表达了对职业生涯的愿望和希望,认为人类有能力不断改善自己。 ... [详细]
  • 如何用UE4制作2D游戏文档——计算篇
    篇首语:本文由编程笔记#小编为大家整理,主要介绍了如何用UE4制作2D游戏文档——计算篇相关的知识,希望对你有一定的参考价值。 ... [详细]
  • XML介绍与使用的概述及标签规则
    本文介绍了XML的基本概念和用途,包括XML的可扩展性和标签的自定义特性。同时还详细解释了XML标签的规则,包括标签的尖括号和合法标识符的组成,标签必须成对出现的原则以及特殊标签的使用方法。通过本文的阅读,读者可以对XML的基本知识有一个全面的了解。 ... [详细]
  • 不同优化算法的比较分析及实验验证
    本文介绍了神经网络优化中常用的优化方法,包括学习率调整和梯度估计修正,并通过实验验证了不同优化算法的效果。实验结果表明,Adam算法在综合考虑学习率调整和梯度估计修正方面表现较好。该研究对于优化神经网络的训练过程具有指导意义。 ... [详细]
author-avatar
吴雨醒
这个家伙很懒,什么也没留下!
PHP1.CN | 中国最专业的PHP中文社区 | DevBox开发工具箱 | json解析格式化 |PHP资讯 | PHP教程 | 数据库技术 | 服务器技术 | 前端开发技术 | PHP框架 | 开发工具 | 在线工具
Copyright © 1998 - 2020 PHP1.CN. All Rights Reserved | 京公网安备 11010802041100号 | 京ICP备19059560号-4 | PHP1.CN 第一PHP社区 版权所有