写个介绍时间序列的，专注于用Python实现时间序列的一些模型及性质，如果想更多得学习理论，推荐参阅一些时间序列的书籍。
参考文献1.《金融时间序列分析》 第2版 Ruey S.Tsay著 王辉、潘家柱 译2.Time Series analysis tsa http://nipy.bic.berkeley.edu/nightly/statsmodels/doc/html/tsa.html
一、与时间序列分析相关的部分基础知识/概念1.1 什么是时间序列简而言之：对某一个或者一组变量x(t)进行观察测量，将在一系列时刻t1,t2,⋯,tn所得到的离散数字组成的序列集合，称之为时间序列。
例如: 某股票A从2015年6月1日到2016年6月1日之间各个交易日的收盘价，可以构成一个时间序列；某地每天的最高气温可以构成一个时间序列。
一些特征:趋势：是时间序列在长时期内呈现出来的持续向上或持续向下的变动。
季节变动：是时间序列在一年内重复出现的周期性波动。它是诸如气候条件、生产条件、节假日或人们的风俗习惯等各种因素影响的结果。
循环波动：是时间序列呈现出得非固定长度的周期性变动。循环波动的周期可能会持续一段时间，但与趋势不同，它不是朝着单一方向的持续变动，而是涨落相同的交替波动。
不规则波动：是时间序列中除去趋势、季节变动和周期波动之后的随机波动。不规则波动通常总是夹杂在时间序列中，致使时间序列产生一种波浪形或震荡式的变动。只含有随机波动的序列也称为平稳序列。
1.2平稳性在百度词条中是这样粗略的讲的：平稳时间序列粗略地讲，一个时间序列，如果均值没有系统的变化（无趋势）、方差没有系统变化，且严格消除了周期性变化，就称之是平稳的。
我们不妨先来看下：
IndexData = DataAPI.MktIdxdGet(indexID=u"",ticker=u"000001",beginDate=u"20130101",endDate=u"20140801",field=u"tradeDate,closeIndex,CHGPct",pandas="1")
IndexData = IndexData.set_index(IndexData['tradeDate'])
IndexData['colseIndexDiff_1'] = IndexData['closeIndex'].diff(1) # 1阶差分处理
IndexData['closeIndexDiff_2'] = IndexData['colseIndexDiff_1'].diff(1) # 2阶差分处理
IndexData.plot(subplots=True,figsize=(18,12))
细心的朋友已经发现，下面两张图，实际上是对第一个序列做了差分处理，方差和均值基本平稳，成为了平稳时间序列，后面我们会谈到这种处理。
下面可以给出平稳性的定义了：
严平稳：
如果对所有的时刻t，任意正整数k和任意k个正整数(t1,t2…tk),
的联合分布与
的联合分布相同，我们称时间序列{rt}是严平稳的。
也就是，
的联合分布在时间的平移变换下保持不变，这是个很强的条件。而我们经常假定的是平稳性的一个较弱的方式
弱平稳：
若时间序列{rt}满足下面两个条件:
则时间序列{rt}为弱平稳的。即该序列的均值，rt与rt−l的协方差不随时间而改变，l为任意整数。
在金融数据中，通常我们所说的平稳序列，是弱平稳的。
差分
回头我们再谈之前说的差分操作：
差分（这里为前向），就是求时间序列{rt}在t时刻的值rt与t−1时刻的值rt−1的差不妨记做dt,则我们得到了一个新序列{dt}，为一阶差分，对新序列{dt}再做同样的操作，则为二阶差分。
通常非平稳序列可以经过d次差分，处理成弱平稳或者近似弱平稳时间序列。回头看上图，我们发现二阶差分得到的序列比一阶差分效果更好。
1.3 相关系数和自相关函数1.3.1 相关系数对于两个向量，我们希望定义它们是不是相关。一个很自然的想法，用向量与向量的夹角来作为距离的定义，夹角小，就距离小，夹角大，就距离大。
早在中学数学中，我们就经常使用余弦公式来计算角度：
我们发现，相关系数实际上就是计算了向量空间中两个向量的夹角！协方差是去均值后两个向量的内积！
如果两个向量平行，相关系数等于1或者-1，同向的时候是1，反向的时候就是-1。如果两个向量垂直，则夹角的余弦就等于0，说明二者不相关。两个向量夹角越小，相关系数绝对值越接近1，相关性越高。 只不过这里计算的时候对向量做了去均值处理，即中心化操作。而不是直接用向量X,Y计算。
对于减去均值操作，并不影响角度计算，是一种“平移”效果，如下图所示：
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
a = pd.Series([9,8,7,5,4,2])
b = a - a.mean() # 去均值
plt.figure(figsize=(10,4))
a.plot(label='a')
b.plot(label='mean removed a')
plt.legend()
相关系数度量了两个向量的线性相关性，而在平稳时间序列{rt}中，我们有时候很想知道，rt与它的过去值rt−i的线性相关性。这时候我们把相关系数的概念推广到自相关系数。
下面给出示例：
from scipy import stats
import statsmodels.api as sm # 统计相关的库
data = IndexData['closeIndex'] # 上证指数
m = 10 # 我们检验10个自相关系数
acf,q,p = sm.tsa.acf(data,nlags=m,qstat=True) ## 计算自相关系数 及p-value
out = np.c_[range(1,11), acf[1:], q, p]
output=pd.DataFrame(out, columns=['lag', "AC", "Q", "P-value"])
output = output.set_index('lag')
output
我们取显著性水平为0.05，可以看出，所有的p-value都小于0.05；则我们拒绝原假设H0。
因此，我们认为该序列，即上证指数序列，是序列相关的
我们再来看看同期上证指数的日收益率序列：
data2 = IndexData['CHGPct'] # 上证指数日涨跌
m = 10 # 我们检验10个自相关系数
acf,q,p = sm.tsa.acf(data2,nlags=m,qstat=True) ## 计算自相关系数 及p-value
out = np.c_[range(1,11), acf[1:], q, p]
output=pd.DataFrame(out, columns=['lag', "AC", "Q", "P-value"])
output = output.set_index('lag')
output
可以看出，p-value均大于显著性水平0.05。我们选择假设H0，即上证指数日收益率序列没有显著的相关性
1.4 白噪声序列和线性时间序列1.4.1 白噪声序列随机变量X(t)(t=1,2,3……)，如果是由一个不相关的随机变量的序列构成的，即对于所有S不等于T，随机变量Xt和Xs的协方差为零，则称其为纯随机过程。
对于一个纯随机过程来说，若其期望和方差均为常数，则称之为白噪声过程。白噪声过程的样本实称成为白噪声序列，简称白噪声。之所以称为白噪声，是因为他和白光的特性类似，白光的光谱在各个频率上有相同的强度，白噪声的谱密度在各个频率上的值相同。
1.4.2 线性时间序列
到目前为止介绍了一些基本知识和概念，如平稳性、相关性、白噪声、线性序列，介绍的过程中并没有太深入，目前来说“够用”了，一些细节会在后面章节补充。下面开始介绍一些线性模型。
二、自回归(AR)模型
下面我们就用该方法，检验上证指数日收益率序列的平稳性
temp = np.array(data2) # 载入收益率序列
model = sm.tsa.AR(temp)
results_AR = model.fit()
plt.figure(figsize=(10,4))
plt.plot(temp,'b',label='CHGPct')
plt.plot(results_AR.fittedvalues, 'r',label='AR model')
plt.legend()
我们可以看看模型有多少阶
print len(results_AR.roots)
17
可以看出，自动生成的AR模型是17阶的。关于价次的讨论在下节内容，我们画出模型的特征根，来检验平稳性
pi,sin,cos = np.pi,np.sin,np.cos
r1 = 1
theta = np.linspace(0,2*pi,360)
x1 = r1*cos(theta)
y1 = r1*sin(theta)
plt.figure(figsize=(6,6))
plt.plot(x1,y1,'k') # 画单位圆
roots = 1/results_AR.roots # 注意，这里results_AR.roots 是计算的特征方程的解，特征根应该取倒数
for i in range(len(roots)):
plt.plot(roots[i].real,roots[i].imag,'.r',markersize=8) #画特征根
plt.show()
可以看出，所有特征根都在单位圆内，则序列为平稳的！
2.2 AR(p)模型的定阶一般有两种方法来决定p：
第一种：利用偏相关函数(Partial Auto Correlation Function,PACF)
第二种：利用信息准则函数
2.2.1 偏相关函数判断p对于偏相关函数的介绍，这里不详细展开，重点介绍一个性质：
AR(p)序列的样本偏相关函数是p步截尾的。
所谓截尾，就是快速收敛应该是快速的降到几乎为0或者在置信区间以内。
具体我们看下面的例子，还是以之前上证指数日收益率序列
fig = plt.figure(figsize=(20,5))
ax1=fig.add_subplot(111)
fig = sm.graphics.tsa.plot_pacf(temp,ax=ax1)
我们看出，按照截尾来看，模型阶次p在110 ，但是之前调用的自动生成AR模型，阶数为17,这一点我也有些不解。。
当然，我们很少会用这么高的阶次。。
2.2.2 信息准则— AIC、BIC、HQ现在有以上这么多可供选择的模型，我们通常采用AIC法则。我们知道：增加自由参数的数目提高了拟合的优良性，AIC鼓励数据拟合的优良性但是尽量避免出现过度拟合(Overfitting)的情况。所以优先考虑的模型应是AIC值最小的那一个。赤池信息准则的方法是寻找可以最好地解释数据但包含最少自由参数的模型。不仅仅包括AIC准则，目前选择模型常用如下准则：
下面我们来测试下3种准则下确定的p，仍然用上证指数日收益率序列。为了减少计算量，我们只计算间隔前10个看看效果。
aicList = []
bicList = []
hqicList = []
for i in range(1,11): #从1阶开始算
order = (i,0) # 这里使用了ARMA模型，order 代表了模型的(p,q)值，我们令q始终为0，就只考虑了AR情况。
tempModel = sm.tsa.ARMA(temp,order).fit()
aicList.append(tempModel.aic)
bicList.append(tempModel.bic)
hqicList.append(tempModel.hqic)
plt.figure(figsize=(15,6))
plt.plot(aicList,'r',label='aic value')
plt.plot(bicList,'b',label='bic value')
plt.plot(hqicList,'k',label='hqic value')
plt.legend(loc=0)
可以看出，3个准则在第一点均取到最小值，也就是说，p的最佳取值应该在1，我们只计算了前10个，结果未必正确。
Whatever，我们的目的是了解方法。
当然，利用上面的方法逐个计算是很耗时间的，实际上，有函数可以直接按照准则计算出合适的order，这个是针对ARMA模型的，我们后续再讨论。
2.3 模型的检验根据式2.0，如果模型是充分的，其残差序列应该是白噪声，根据我们第一章里介绍的混成检验，可以用来检验残差与白噪声的接近程度。
我们先求出残差序列：
delta = results_AR.fittedvalues - temp[17:] # 残差
plt.figure(figsize=(10,6))
#plt.plot(temp[17:],label='original value')
#plt.plot(results_AR.fittedvalues,label='fitted value')
plt.plot(delta,'r',label=' residual error')
plt.legend(loc=0)
然后我们检查它是不是接近白噪声序列
acf,q,p = sm.tsa.acf(delta,nlags=10,qstat=True) ## 计算自相关系数 及p-value
out = np.c_[range(1,11), acf[1:], q, p]
output=pd.DataFrame(out, columns=['lag', "AC", "Q", "P-value"])
output = output.set_index('lag')
output
观察p-value可知，该序列可以认为没有相关性，近似得可以认为残差序列接近白噪声。
2.4 拟合优度及预测2.4.1 拟合优度我们使用下面的统计量来衡量拟合优度：
但是，对于一个给定的数据集，R2是用参数个数的非降函数，为了克服该缺点，推荐使用调整后的R2:
它的值在0-1之间，越接近1，拟合效果越好。
下面我们计算之前对上证指数日收益率的AR模型的拟合优度。
score = 1 - delta.var()/temp[17:].var()
print score
0.0405231650285
可以看出，模型的拟合程度并不好，当然，这并不重要，也许是这个序列并不适合用AR模型拟合。
2.4.2 预测我们首先得把原来的样本分为训练集和测试集，再来看预测效果，还是以之前的数据为例：
train = temp[:-10]
test = temp[-10:]
output = sm.tsa.AR(train).fit()
output.predict()
predicts = output.predict(355, 364, dynamic=True)
print len(predicts)
comp = pd.DataFrame()
comp['original'] = temp[-10:]
comp['predict'] = predicts
comp
10
